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I.  INTRODUCTION 

Two  methods  for  calculating  electromagnetic  plane  wave  scattering 
from  a  perfectly  conducting  TM  cylindrical  surface  with  a  narrow  but 
infinitely  long  slot  were  described  in  [ 1 ] .  These  two  methods  are  the 
method  of  solution  with  pseudo-image  and  the  method  of  solution  without 
pseudo- image .  A  Fourier  series  method  of  solution  for  calculating  the 
scattering  of  an  electromagnetic  plane  wave  normally  incident  on  a  narrow 
but  infinitely  long  slot  in  a  perfectly  conducting  TM  circular  cylindrical 
surface  was  also  described  in  [1].  In  what  follows,  two  computer  programs 
are  described  and  listed,  one  that  implements  the  methods  of  solution  with 
and  without  pseudo-image,  and  one  that  implements  the  Fourier  series  method 
of  solution. 

The  computer  program  for  the  methods  of  solution  with  and  without 
pseudo- image  consists  of  a  main  program  and  the  subroutines  BES,  DECOMP, 
and  SOLVE.  The  subroutine  BES  is  described  and  listed  in  Section  II,  the 
subroutines  DECOMP  and  SOLVE  in  Section  III,  and  the  main  program  in 
Section  IV.  The  main  program  calculates  the  magnetic  current  coefficient 
V  of  [1,  eq.  (9)]  by  the  methods  of  solution  with  and  without  pseudo-image. 

The  method  of  solution  with  pseudo-image  is  described  in  [1,  Sections  III 
to  VII].  The  method  of  solution  without  pseudo-image  is  described  in  [1, 
Appendix  A].  The  main  program  reads  input  data  from  the  file  MAUTZ1.DAT 
and  writes  output  data  on  the  file  MAUTZ2.DAT.  In  Section  IV,  the  contents 
of  the  files  MAUTZ1.DAT  and  MAUTZ2.DAT  are  listed  when  the  input  and  output 
data  are  for  the  example  of  [1,  Section  VIII]. 

The  computer  program  for  the  Fourier  series  method  of  solution  consists 
of  a  main  program  and  the  subroutines  BES,  BESJ1 ,  and  BESJY .  The  subroutine 
BES  is  the  same  as  the  one  mentioned  in  the  previous  paragraph.  The  sub¬ 
routine  BESJ1  is  described  and  listed  in  Section  V,  the  subroutine  BESJY 


in  Section  VI,  and  the  main  program  in  Section  VII.  The  main  program 
obtains  an  approximation  for  the  magnetic  current  coefficient  V  of  [1, 
eq.  (B— 25) ]  by  truncating  the  infinite  series  in  [1,  eq.  (B— 25) ]  at 
n  ■  N2  where  N2  is  entered  as  input  data.  The  main  program  reads  input 
data  from  the  file  MAUTZ3.DAT  and  writes  output  data  on  the  file  MAUTZ4.DAT. 

In  Section  VII,  the  contents  of  the  files  MAUTZ3.DAT  and  MAUTZ4.DAT  are 
listed  when  the  input  and  output  data  are  for  the  example  of  [1,  Section  VIII] 
and  N2  =  10,000. 


II.  THE  SUBROUTINE  BES 


The  subroutine  BES(N,X)  puts 


J  (X)  in  B J (n+1 )~) 


,  n  ■  0,  1,2, ...  ,N 


(1) 


Yn(X)  in  BY (n+l)J 


where  J  is  the  Bessel  function  of  the  first  kind  of  order  n,  and  Y  is 
n  n 


the  Bessel  function  of  the  second  kind  of  order  n.  Moreover,  N  is  a  non¬ 


negative  integer  and  X  is  a  non-negative  real  number.  The  subroutine  BES 
is  listed  at  the  end  of  this  section.  In  the  COMMON  statement  on  line  3 


of  this  listing,  YO,  PI2,  PI4,  and  PI7  are  input  variables  and  BJ  and  BY 
are  the  output  variables  that  appear  in  (1).  The  elements  of  YO  will  be 
described  later.  Their  values  should  not  deeply  concern  the  user  because 


he  will  never  have  to  change  them.  The  values  of  the  input  variables 
Y0(1),  Y0 (2) .... ,Y0 (33)  are  listed  as 


0.  3072582E*04  0. 736S758E*04-C. 6C85 1 00E*03  0.17102347*02-0. 2 27 10017* 00 
0.  1600171  E-02-C.5961089E-05  0. 9  545  7  73  E-08  0. 4 1 6 3 1 50E* 0 5  C .  3  4202  1  IE*  C  "> 
0.  1C00C CO E*0 1-0. 602472 7E *04  0.1613512E*04-C.7532210E*02  0.  14025  9CE  +  01 
0.  12756C2E-01  0 . 583 2787E-0 4-0.  1  1 07 6 99 E-06  0. 30729  46E*05  0 .  28 864  3  1 E* C3 
0.  10C00C0E  +  0  1  C.99S9599E*00-C.1097659E-02  0. 246 1  455E-04  0 .  1 C0000 C E  +  0  1 
0.  1829893 E- 02-0. 319  1328 E-0 4-0.  1562498E-0  1  0.  14 27 0 79E-0 3- 0 . 5 9 374 34E- 05 
0.4667498E-01-0.1998720E-03  C. 7317495 E- 05 


where  the  jth  number  in  the  ith  row  is  YO(5(i-l)+j)  for  j=l,2,.... 


5  and 


i=l , 2 , 3 , .  .  .  ,  7 .  The  values  of  the  remaining  input  variables  in  the  common 
statement  are  given  by 

PI2  =  2h  ! 

i 

V  PI4  =  m/4  >  (2) 

I 

PI7  =  3-/4  \ 

Minimum  allocations  are  given  by 

COMMON  BJ (N+l) ,  BY (N+l) 

DIMENSION  AJ (ll+N+2 [X] ) 

where  [X]  is  the  largest  integer  that  does  not  exceed  X. 

The  Bessel  functions  {j  (X),  n=0, 1 , 2 , . . . ,N  ‘  are  calculated  in 

n 

lines  5  to  23.  As  suggested  in  [2,  Sec.  9.12.,  Example  1],  J^(X)  and 

-20 

Jmz_iW  are  set  equal  to  the  arbitrary  values  of  zero  and  10  ,  respec¬ 

tively,  where 

MZ  =  10  +  N  +  2 [X]  (3) 

This  is  done  in  lines  5  to  8  where  .1  and  J^^fX)  are  stored  in 

AJ (MZ+1)  and  AJ(MZ),  respectively.  In  DO  loop  16,  the  recurrence 
relation  [2,  eq .  (9.1.27)] 

jn(x)  .IMUjn^CX)  -  J„+2(X,  (4) 

is  used  to  calculate  (j  (X),  n  =  MZ-2,  MZ-3,...,0K  Line  13  puts 

n 

^MK  L  ^  in  AJ  (MK)  .  According  to  [2,  eq  .(  9 . 1 .46)  ],  the  calculated 

values  of  {J  (X)]  have  to  be  normalized  bv  dividinc  bv  >  where 
n 

x  =  J  (X)  +  2J~ (X)  +  2J. (X)  +  ...  (5) 

o  2  4 

Each  of  the  calculated  values  of  J  (X),  n=0,l,2,  ...,N  is  divided  by  t  in 


DO  loop  17. 


’viar.  v  "c  *  -  xa 


The  Bessel  functions  Y  (X),  n=0, 1 , 2 , . . . , X '  are  calculated  in 

n 

lines  24  to  68.  For  X  <  8,  we  let  [3,  eq.  (6.8.27)] 


where 


Y  (X)  =  Y  (X)  +  -  J  (X)  In  (X) 
O  O  -T  o 


l  (POm)X 

Y  (X)  =  ^ - 

o  2 


l  (QOm)X 
m=0 


where  {POm}  and  {QOm}  are  the  coefficients  given  under  the  heading 
YZERO  6234  on  page  309  of  [3].  For  X  <  8,  we  also  let  [3,  eq.  (6.8.28)] 


YX(X)  =  Y:(X)  +  -  fJL(X)  ln(X)  -  1/X] 


where 


Y  (v)  l  (P0m)XZm 

=  5=2 _  (9) 

X  2  K  ’ 

l  (Q0m)X2m 
m=0 

where  ■ POm)  and  ;QOm}  are  the  coefficients  given  under  the  heading  YONE 

6434  on  page  315  of  [3].  The  coefficients  { POm I  of  (7),  \Q0tn)  of  (7),  ^POm 

of  (9),  and  'QOm)  of  (9)  reside  in  YO(l)  to  Y0(8),  Y0(9)  to  Y0(11),  Y0(12) 

to  Y0(18),  and  Y0(19)  to  Y0(21),  respectively.  Lines  41  and  42  put  Yq(X) 

of  (7)  in  S,  and  line  44  puts  Y  (X)  of  (6)  in  BY(1).  Lines  46  and  47 

o 

put  Y  (X)/X  of  (9)  in  S,  and  line  48  puts  Y^X)  of  (8)  in  BY(2). 

For  X  '•  8,  we  let  [3,  eq.  (6.8.17)] 


Y  (X)  = 
n 


~  ■?  (X)  sin (X  )+0  (X)  cos (X  );,  n=0,l 

X  n  n  n  n 


In  (10), 


X 


and 


where 


X  = 
n 


(2n+l)  . 

4 


P  (X)  =  l  (POm)Z2m 
0  m=0 


Z  =  8/X 


(ID 


(12) 


(13) 


and  ; POm,  are  the  coefficients  given  under  the  heading  PZERO  6501  on 
page  320  of  [3],  Moreover, 

2  7 

Q  (X)  =  Z  )  (POm)Z  (14) 

o  u 

m=0 

where  {POm}  are  the  coefficients  given  under  the  heading  QZERO  6900  on 
page  327  of  [3],  Regarding  (14),  there  are  errors  on  page  149  of  [3]  . 

On  that  page,  each  Qq(x)/x  should  be  replaced  by  Q^ix) / z  and  each  Q^(x)/x 
by  Q1(x)/z.  In  (10) , 

2  7 

P  (X)  =  l  (POm)Z  (15) 

m=0 

where  {POm-  are  the  coefficients  given  under  the  heading  PONE  6701  on 
page  323  of  [ 3 ] ,  and 

2 

Q  (X)  =  Z  T  (POm) ZZm  (16) 

m=0 

where  ‘POm}  are  the  coefficients  given  under  the  heading  QON'E  /101  on  page 

331  of  [3].  The  coefficients  •'  POm  •  of  (12),  POm’  of  (15),  ’POm’  of  (14), 

and  ^POrnf  of  (16)  reside  in  YO(22)  to  YO(24),  YO(25)  to  YO(27),  YO(28) 

to  Y0(30) ,  and  Y0(31)  toYO(33),  respectively.  Line  54  puts  P  (X)  of  (12) 

in  SI,  line  55  puts  Q  (X)  of  (14)  in  S3,  line  56  puts  X^  of  (11)  in  S5, 

and  line  57  puts  Y  (X)  of  (10)  in  BY(1).  Line  59  puts  P  (X)  of  (151  in 
o  t 

S2,  line  60  puts  (X)  of  (16)  in  S4,  line  61  puts  X^  oi  (11)  in  S6,  and 
line  62  puts  Y.  (X)  of  (10)  in  BY(2). 


If  N  >  1,  then  DO  loop  21  uses  the  recurrence  relation  [2,  eq 


(9.1.27)] 


Y  (X)  =  2  Y  (X)  -  Y  .(X) 

n  X  n-1  n-2 

to  put  Y  (X)  in  BY(n+l)  for  n»2,3,...,N. 
n 


C01  c 


005 
C06 
C07 
C  08 
009 
CIO 
Cl  1 
C 1  2 
013 
C 1  4  16 
C  1  5 
Cl  6 
C 1 7 
C 1  8  15 
C 1  9 
C2C 
C21 
022 
C  2  3  17 


LISTING  CP  THE  SOBBCUTINE  BES 
SOBBOtJTINE  EES  ( N ,  X) 

COMMON  10(33)  ,PI2,PI4,PI7,BJ  (100)  , BY  (10  0) 
DIMENSION  A  J  (100) 

NZ=1C*N*2*IFIZ(X) 

IF  (X.  IT.  1.  E-3)  MZ=4+N 
A  J  (NZ  ♦  1 )  =0. 

AJ (HZ)=1.E-20 
M1  =  MZ-  1 
X2=2./X 
DO  16  K=  1,M  1 
BK=PZ-K 

AJ (MK)=X2*FLCAT  (KK)  *AJ(MK+1) -AJ  (MS*?) 

CONTINUE 

ALP=«  5*A J  (1) 

EC  15  J=3,MZ,2 
ALP=ALP*AJ  (J) 

CCNTINDE 

ALP=2-*AL? 

NE=N+1 

DO  17  K=  1 »  NP 
EJ  <K)=AJ  (K)  /ALP 
CONTINOE 


C2  4 
C25  18 
C26 
C27 
C2d 
C29 
C3Q 
031 
C32 
033 
034 
C35 
C36 
037 
C38 
C39 
C4C 
041  22 
C42 
C4J 
C44 
045 
C46 
047 
C4  8 
049 
CSO  19 
051 
C  52 
053 
054 
055 
056 
057 
058 
C59 
060 
06  1 
062 
063  20 
064 
065 
066 
067  21 
068 
C69 


1.  E-6) 


TC  22 


TC  22 


E-4)  GC  TC  22 


IF(X-8.)  18,18,19 

Z  1=X*X 

22=0. 

Z  2=0. 

24=0. 

Z  5=0. 

Z6=C. 

Z7=  0. 

IF  (Zl.  LT.  1.  E-1  1)  GO  TC  22 
Z2=Z  1*Z 1 
Z3=Z2*Z  1 

IF  (Z  1.  IT.  1.  E-6)  GC  TC  22 
24  =  Z  3*2 1 
Z  5  =  Z  4  *  Z 1 

IF  (Z  1.  IT.  1.  E-4)  GC  TC  22 
Z6=Z5*Z 1 
Z7=Z6*Z 1 

S=  (YC  (1)  ♦  YO  (2)  *Z1*Y0  (3)  *22  + YO  (4)  *Z3*Y0  (5)  *Z4  +  Y0  (6)  *Z5* 
YO  (7)  *Z6 ♦ YO  (8)  *Z7)  /  (YO  (9)  *Y0  ( 10)  *ZHY0  ( 1  1)  *Z2) 

TLG  =  AICG  ( X) 

EY  <1)=S*PI2*3J  (  1)  *TIG 
IF(N-IE.C)  PETOSN 

S=(YO  (12)  >Y0(13)  *ZUY0(14)  *Z2  +  Y0(15)  *Z3  +  Y0  I  1 6) *Z4  + 

Y0(  17)  *Z5*Y0  (18)  *26)  /  (YO  (19)  +YO  (20)  *Z1  +  Y3  (21)  *Z2) 

EY  (2)  =X*S*PI2*  (BJ  (2)  *TIG-  1./X) 

GO  TC  20 
2=9. /X 
Z1=Z*Z 
Z2  =  Z  1*Z  1 
S  =  SQRT  (PI2/X) 

S1=  YO  (22)  *Y3  (2  3)  *Z  1  +  YO  (24)  *Z2 
S3  =  Z*(YC(28)  *■  Y  0  (29)  *Z1*Y0(3G)  *Z2) 

S  5  =  X- FI  4 

EY(1)=S*(S1*SIN(S5)*S3*CCS(£5)) 

IF(N.IE.O)  HETTJBN 

S2  =  Y0  (2  5)  ♦  Y  0  (26)  *Z1  +  Y0{2  7)  *Z2 

S4=Z*  (YO  (31)  *Y0  (3  2)  *Z1*Y0  (3  3)  *Z2) 

S6=  X- F 17 

EY  (2)  =  S*  (S2  *SI N  (S6)  ♦ S4*CCS  (S6)  ) 

IF  (NE.  IE.  2)  FETUPN 
DO  21  K=  3, NP 
F2=  K-2 

BY  (K)  =  X2*FICAT  (K2)*EY(K-1)-EY(K2) 

CONTINUE 

RETURN 

END 


III.  THE  SUBROUTINES  DECOMP  AND  SOLVE 

The  subroutines  DECOMP (N,  IPS,  UL)  and  SOLVE (X,  IPS,  UL,  B,  X)  solve 
a  system  of  linear  equations  in  N  unknowns.  The  input  to  DECOMP  consists 
of  N  and  the  N  by  N  matrix  of  coefficients  on  the  left-hand  side  of  the 
matrix  equation  stored  by  columns  in  UL.  The  output  from  DECOMP  is  IPS 
and  UL.  This  output  is  fed  into  SOLVE.  The  rest  of  the  input  to  SOLVE 
consists  of  N  and  the  column  of  coefficients  on  the  right-hand  side  of  the 
matrix  equation  stored  in  B.  SOLVE  puts  the  solution  to  the  matrix  equa¬ 
tion  in  X. 

Minimum  allocations  are  given  by 
COMPLEX  UL(N*N) 

DIMENSION  SCL(N),  IPS(N) 

in  DECOMP  and  by 

COMPLEX  UL(N*N) ,  B(N) ,  X(N) 

DIMENSION  IPS (N) 

in  SOLVE. 

More  detail  concerning  DECOMP  and  SOLVE  is  on  pages  46-49  of  [41 . 


CO  1  c 

LISTING  OF  THE  SOBBCOTINE  DECCMP 

002 

SUBROUTINE  DECOMP ( N, IPS, UL) 

C03 

C  CM  PIS I  UL(loOO)  , PIVOT, EK 

004 

DIMENSION  SCL  (40)  ,  IPS  (40) 

005 

DC  5  1=1, N 

006 

IPS  (I)  =1 

007 

RN=0. 

008 

J  1=1 

009 

DO  2  J=  1  ,  N 

010 

OLH=AES  (REAL  (01  ( J1)  )  )  +  ABS  ( AIMAG  (UL  (J1)  )  ) 

Oil 

J 1= J  1  +  N 

012 

IF  (RN-ULM)  1,2,2 

013  1 

RN=  OLfl 

014  2 

CONTINOE 

01  5 

SCL  (I)  =  1./RN 

016  5 

CONTINOE 

017 

NH1=N-1 

018 

K2=0 

019 

DO  17  K=1 ,N  Ml 

020 

EIG=  0. 

021 

DO  11  I=K,N 

022 

IP=IPS  (I) 

023 

IPK=IP+K2 

02  4 

SIZE=(ABS  (REAL  (OL(IPK) )  )  ♦  ABS  (AIM AG  (UL  (IPK)  )  )  )  *SCL  (IP) 

025 

IF  (SIZE-BIG)  11,11,10 

026  10 

EIG=SIZE 

027 

IPV=I 

028  11 

CONTINOE 

02  9 

IF(IEV-K)  14,15,14 

030  14 

J=IPS(K) 

031 

IPS  (K)  =IPS  (IPV) 

032 

IPS(IPV)  =J 

033  15 

K  PP  =  I PS (K) ♦ K2 

034 

PI VOT  =  OL  (KPP) 

C3  5 

KP 1  =  K>  1 

036 

CO  16  I=KP  1 , N 

03  7 

KP=KFP 

038 

IP=IPS  (I)  +K2 

C39 

EH=-OL(IF)/PIVCT 

040  18 

OL (IP) --EH 

041 

DC  16  J=KP1,N 

042 

IP=IP+N 

043 

KP=KF*N 

044 

OL  ( IP)  =OL  (IP)  ♦  EM*OL(KP) 

C45  16 

CCNTINOE 

04  b 

K2=K2*N 

04  7  17 

CCNTINOE 

048 

RETORN 

C49 

ENO 

LISTING  OF  THE  SUBROUTINE  SOLVE 
SUBROUTINE  SOLVE (N, I?S,UL, 3,X) 
COMPLEX  UL  {  1600)  ,  E  (  4C)  ,  X(4  0)  ,SUM 
0I3ZNSICM  IPS  (40) 

N  P  1=  N  1 
IP=IPS  (1) 

X  (1)  —  B  (IP) 

DO  2  1  =  2  ,  N 
IP=IPS (I) 

I P  B=  I P 
Ifl  1=1-1 
SOM=0. 

CO  1  J=  1,13  1 
SUM=S03+0L(IF)*X  (J) 

I?=I?+N 

X  (I)=B  (I PB)  -SOW 
K  2=N*  (N-  1) 

IP=IES  (N) ♦K2 
X  (N)  =  X  { N)  /OL  {I  P) 

DO  4  IB  AC  R=  2  ,  H 
I  =  NP 1-IEACK 
S2=K2-N 
IPI=IPS  (I)  *K2 
IP1=I+1 
503=0. 

IP=IPI 

CO  3  J=IP1,N 
IP=IP+N 

SUH  =  SOM*OL(IP)  *X(J) 

X  (I)  =  (X  (I)-SOH)/OL  (IPI) 

RETURN 


IV.  THE  MAIN  PROGRAM  FOR  THE  METHODS  OF  SOLUTION  WITH  AND  WITHOUT 
PSEUDO-IMAGE 

The  main  program  for  the  methods  of  solution  with  and  without 
pseudo-image  calculates  the  magnetic  current  coefficient  V  of  [1,  eq. 

(9)]  by  means  of  the  methods  of  solution  with  and  without  pseudo- image. 
The  method  of  solution  with  pseudo-image  is  described  in  [1,  Sections 
III  to  VII] .  The  method  of  solution  without  pseudo-image  is  described 
in  [1,  Appendix  A].  Input  data  are  read  from  the  file  MAUTZ1.DAT,  output 
data  are  written  on  the  file  MAUTZ2.DAT,  and  the  subroutines  BES,  DECOMP, 
and  SOLVE  are  called.  The  subroutine  BES  was  described  and  listed  in 
Section  II  and  the  subroutines  DECOMP  and  SOLVE  in  Section  III. 

The  main  program  for  the  methods  of  solution  with  and  without 
pseudo-image  is  listed  at  the  end  of  this  section.  In  the  open  state¬ 
ments  on  lines  9  and  10  of  this  listing,  MAUTZ1.DAT  is  the  input  data 
file  and  MAUTZ2.DAT  is  the  output  data  file.  The  input  data  are  read 
early  in  the  program  according  to 

READ (20,  11)  N,  NG,  BK,  PINC 
11  FORMAT (213 ,  2E14.7) 

READ (20,  13) (X(I) ,  I  =  1 ,N) 

13  FORMAT (5E14. 7) 

READ (20 ,  13) (Y (I)  ,  I  =  1,N) 

READ (20,  13) (YO(I) ,  1=1,  33) 

READ(20 ,  13) (XG(I) ,  I  =  1,  NG) 

READ(20,  13) (AG (I) ,  I  =  1 ,  NG) 


Here,  N  is  the  integer  that  appears  in  [1,  eq .  (22)].  That  is,  N  is  the 
number  of  electric  current  expansion  functions  on  the  complete  conducting 


surface  S  .  An  NG-point  Gaussian  quadrature  formula  is  used  whenever 
Gaussian  quadrature  is  called  for.  Specifically,  the  integral  from  -1 
to  1  with  respect  to  x  of  a  function  f (x)  is  approximated  by  [5, 
Appendix  A] 


NG 


f(x)dx  =  l  A^NGWNG)) 


(18) 


-1 


1= 1 


where  X^G^  and  A^G^  are,  respectively,  the  abscissa  and  weight  given 


in  [5,  Appendix  A].  Still  in  the  first  read  statement,  BK  is  the  wave 


number  k  that  appears  in  [1,  eq.  (44)],  and  PINC  is  the  incident  angle 
inc 


in  [1,  eq.  (44)].  In  the  next  two  read  statements,  (X(I) ,  Y(I))  are 
the  (x,y)  coordinates  of  the  point  P  in  [1,  Fig.  1],  That  is,  X(I)  and 


Y(I)  are,  respectively,  the  quantities  x^  and  y^.  in  [1,  eq .  (53)].  Here, 


BK  is  in  reciprocal  meters,  PINC  is  in  radians,  and  (X(I) ,  Y(I))  are  in 
meters.  The  array  YO  contains  input  data  for  the  subroutine  BES.  The 
values  of  the  elements  of  YO  were  given  in  Section  II.  These  values 
should  not  deeply  concern  the  user  because  he  will  never  have  to  change 


them.  In  the  last  two  read  statements,  XG(I)  and  AG(I)  are,  respectively. 


the  abscissa  X^G^  and  weight  A^G^  in  (18). 


Minimum  allocations  are  given  by 


COMMON  BJ (2) ,  BY ( 2 ) 

COMPLEX  Z (N*N) ,  VA(N) ,  VINC(N),  CURI(N),  CURA(N) 

DIMENSION  X (N) ,  Y (N) ,  XG(NG) ,  AG(NG) ,  DX(N) 

DIMENSION  DY (N) ,  G(N),  XP(N),  YP(N),  IPS(N) 

In  the  present  program,  BJ(2)  and  BY(2)  are  sufficient.  However,  BJ(100) 
and  BY (100)  were  used  in  the  listing  of  this  program.  Tin’  reason  for  this 


is  explained  in  the  next  two  sentences. 


The  allocations  for  BJ  and  BY  in 


the  present  program  must  be  the  same  as  the  allocations  in  the  subroutine 

BES.  The  same  subroutine  BES  was  also  used  with  the  program  of  Section  VII, 

and  during  that  usage  BJ(2)  and  BY(2)  did  not  suffice. 

Immediately  after  the  main  program  at  the  end  of  this  section,  the 

contents  of  the  input  data  file  MAUTZ1.DAT  and  the  output  data  file 

MAUTZ2.DAT  are  listed  when  the  input  and  output  data  are  for  the  example 

of  [1,  Section  VIII].  The  output  data  file  MAUTZ2  .DAT  contains  all  the 

data  put  out  by  the  write  statements  in  the  main  program  for  the  methods 

of  solution  with  and  without  pseudo-image.  The  contents  of  the  output 

data  file  MAUTZ2.DAT  are  described  in  the  next  two  paragraphs. 

The  input  data  are  written  out  immediately  after  they  are  read  in. 

The  four  numbers  in  the  ith  line  under  the  heading  "VINC"  are,  from  left 

to  right,  ReCAkV^^),  Im^kV*”^) ,  Re (4kVaaC) ,  and  Im(4kVa,?c)  where 
inc  inc 

V„ .  ,  and  V_  .  are  the  (2i-l)th  and  (2i)th  elements  of  the  column  vector 
2i-l  2i 

Vinc  that  appears  in  [1,  eq.  (23)].  Moreover,  "Re"  denotes  real  part, 

"Im"  denotes  imaginary  part,  and  k  is  the  wave  number  that  appears  in 

[1,  eq.  (44)].  Similarly,  the  numbers  written  under  the  heading  "VA"  are 

the  elements  of  4kV3,  those  written  under  the  heading  "Z"  are  the  ele- 
4k 

ments  of  —  Z  , ,  those  written  under  the  heading  "CURI"  are  the  elements 

r  ,  1 

of  nIlnC,  and  those  written  under  the  heading  "CURA"  are  the  elements  of 
Dla.  Here,  Va  appears  in  [1,  eq.  (34)],  Z  denotes  the  first  column  of 

9  -L 

the  square  matrix  Z  that  appears  in  [1,  eq.  (23)],  ”  is  the  intrinsic 
impedance  that  appears  in  [1,  eq.  (45)],  iinc  appears  in  [1,  eq .  (23)1, 
and  fJ  appears  in  (1,  eq.  (34)]. 

Other  output  variables  are  called  CA,  CP,  CIN'C,  YHS,  YAB,  TI,  V, 

YABW,  and  VW.  They  are  conspicuously  identified  in  the  output  data  file 
MAUTZ2.DAT.  For  instance,  the  real  and  imaginary  parts  of  CA  are  preceded 


1 


by  "CA=".  The  output  variable  CINC,  being  real,  is  a  single  number 
preceded  by  "CINC=".  Each  output  variable  mentioned  in  the  first  sentence 
of  this  paragraph  corresponds  to  a  quantity  in  [1].  This  correspondence 
is  revealed  in  Table  1. 


Table  1.  Computer  program  output  variables  and  their  corresponding 
quantities  in  [1]  or  in  the  text. 


Output 

variable 


Corresponding  quantity 
in  [1]  or  in  the  text 


Equation(s)  where  the 
uantity  appears 


N 

N 

[1, 

eq. 

(22)] 

NG 

NG 

(18) 

BK 

k 

[1, 

eq. 

(44)] 

PINC 

^inc 

r 

[1, 

eq. 

(44)] 

X(I) 

xr 

[1, 

eq. 

(53)] 

Y  (I) 

yi 

[1, 

eq. 

(53)] 

YO 

POm  and  QOm 

(7), 

(9), 

(12), (15), (14), (16) 

XG(I) 

y(ng) 

I 

(18) 

AG  (I) 

(NG) 

AI 

(18) 

YINC(I) 

4kV^nC 

[1, 

eq. 

(23)] 

YA  ( I ) 

4  k  Vj 

U, 

eq  . 

(34)] 

Z  ( I ) 

it  z 
•  11 

[1, 

eq. 

(23)] 

CA 

CA 

fl. 

eq  . 

(86)1 

CP 

c. 

[1. 

eq. 

(86)] 

CINC 

krClnC 

[1, 

eq. 

(116)] 

CURI(I) 

inc 

T,II 

[1. 

eq. 

(23)] 

CURA (I ) 

nI? 

[1. 

eq. 

(34)1 

YHS 

2k^Yhs 

[1, 

eq. 

(76)] 

YAB 

4k- (Ya  +  Yb) 

11, 

eq  . 

(10)1 

TI 

4knl 

[1. 

eq. 

(10)] 

V 

V 

[1, 

eq. 

(10)] 

YABW 

4kn(Yaw  +  Ybw) 

fl. 

eq  . 

(A-25) ] 

VW 

vw 

fl. 

eq. 

(A-26 ) 1 

Lines  32  and  33  put  cos  :  and  sin  :  of  fl,  eq .  (44)1  in  CSP 

and  SN’P ,  respectively.  Lines  35  to  37  define  the  common  variables  (2) 
that  are  input  data  for  the  subroutine  BES.  With  regard  to  [1,  eqs. 

(57)  and  (58)],  DO  loop  37  puts  0.5k(xJ+1  -  x^)  in  DX(J)  and  0.5k(y^+j  -  v^) 
in  DY(J).  DO  loop  37  also  puts  0.5)  ^  in  G(J),  kx^  in  XP(J),  and  ky^  in 

YP (J)  where  Yj»  x+,  and  y+  are  given  by  [1,  eqs.  (53),  (55),  and  (56)]. 

4k 

In  nested  DO  loops  38  and  39,  —  is  put  in  Z ( (J- 1) *N+I)  where 

Z  j  is  given  by  [1,  eq.  (52)].  Thus,  as  the  indices  I  and  J  traverse 

4k 

their  ranges,  the  upper  triangular  portion  of  the  square  matrix  —  Z  of 
[1,  eq.  (23)]  is  stored  by  columns  in  the  singly  dimensioned  array  Z. 

When  J-N,  DO  loop  39  puts  4kV^nc  in  VINC(I)  and  4kV^  in  VA(I)  where  V^nC 
and  Vj  are  given  by  [1,  eqs.  (67)  and  (74)]. 

4k 

The  calculation  of  the  off-diagonal  elements  of  —  Z  of  fl,  eq.  (23)] 
is  described  in  this  paragraph  and  the  next  two  paragraphs.  Replacing  ij 
by  IJ  in  [1,  eq.  (52)],  multiplying  both  sides  of  that  equation  by  4k/o  , 
changing  the  variables  of  integration  therein  to  make  the  ranges  of  inte¬ 
gration  from  -1  to  1,  and  then  approximating  both  integrals  by  means  of 
(18),  we  obtain 


4k 


It',,  NG  NG 


I'J  V  .(NG)  v  *(NG)  u(2) 


n  IJ  4 


I  V  >  h^'(*ttot) 


(19) 


K=1 


L=1 


o  v ' IJKL' 


where 


IJKL 


=  /(XX)2  +  (YY)2 


(20) 


where 


XX  =  kxt  -  kxt  +  I  (x T  ,  .  -  xT)xfNG)  -  ^  (x  T , .  -  xT)XT(NG) 


(21) 


Expressions  (21)  and  (22)  were  obtained  by  using  [1,  eqs.  (57)  and  (58)] 

to  express  the  trigonometric  functions  in  [1,  eq.  (54)]. 

Nested  DO  loops  40  and  41  perform  the  double  summation  in  (19).  The 

DO  loop  indices  K  and  L  obtain,  respectively,  the  summation  indices  K  and 

L  in  (19).  DO  loop  41  accumulates  in  Z2  the  summation  with  respect  to  L 

in  (19) .  DO  loop  40  accumulates  in  Z1  the  summation  with  respect  to  K 

in  (19).  Line  80  puts  :>  of  (20)  in  XI.  Depending  on  whether  J  =  N, 

i  J  K-L 

either  line  82  or  line  84  puts  the  Bessel  functions  J  (a  1  and  Y 

o  1JKL  o  IJKL 

in  BJ(1)  and  BY(1),  respectively.  According  to  [2,  eq .  (9.1.4)], 

H(2)(x)  =  J  (x)  -  jY  (x)  (23) 

n  n  n 

Line  86  uses  (23)  with  n=0  and  x=a  T  to  properly  increment  Z2.  Line  93 

I J  PvL 

puts  the  right-hand  side  of  (19)  in  Z( (J-1)*N+I) . 

If,  in  (19),  I  and  J  are  interchanged,  and  the  order  of  summation  is 
interchanged,  the  right-hand  side  of  (19)  remains  unchanged.  Therefore, 

Z  of  (19)  is  svmmetric,  that  is. 


(24) 


Using  (24),  line  102  puts  —  Z  of  (19)  in  Z ( ( T-l ) *N+.T )  .  Thus,  as  I 

and  J  traverse  their  ranges,  line  102  stores  by  columns  the  lower  fri- 

4k 

angular  portion  of  the  square  matrix  —  Z  of  [1,  eq.  (23)]  in  the  singly 
dimensioned  array  Z. 

Multiplying  both  sides  of  [1,  eq .  (67)]  bv  4k,  replacing  i  bv  I 
therein,  expanding  cos  ( 1 ^  -  cinC),  and  using  fl»  eqs.  (57)  and  (58)]  to 


express  the  resulting  cos  i  and  sin  : 


terms,  we  obtain 


t' 


4'l e 


jk(x  cos  ; 


y  sin  : 


(iia-ilii) 

b  x 


where 


SI  -  0.5k  (xI+1  -  Xj.)  cos  4>lnC  +  0.5k(yI+^  -  y^)  sin  ilnC  (26) 

Equation  (25)  is  written  with  the  understanding  that  sin(Sl)/Sl  should  be 
replaced  by  1  whenever  SI  =  0.  If  J  =  N,  then  line  96  inside  nested  DO 
loops  38  and  39  puts  the  right-hand  side  of  (26)  in  SI.  Line  100  puts  the 
right-hand  side  of  (25)  in  VINC(I). 

Since  the  coordinates  (x^,  yN)  and  (xMj_1  ,  yMx1 )  in  [1]  are  given  by 


N+i  'N+r 


(V  V  =  (0,  _w) 


(xN+l’  W  =  (°’  w) 


it  is  evident  from  [1,  eqs.  (55)  and  (56)]  that 


=  0  ; 


kyN  ’  0  J 


It  is  also  evident  from  [1,  eq.  (53)]  that 


k<XN+l  -  V  ‘  0 


k(yN+l  -  V  -  YN 


Taking  advantage  of  (28)  and  (29),  lines  107  to  110  put  4kV  of  (25) 
in  VINC(N). 

It  is  evident  from  (27)  and  (29)  that 


K 


© 


Multiplying  both  sides  of  [1,  eq.  (74)3  by  4k,  replacing  i  by  I  therein, 
changing  the  variables  of  integration  to  make  the  ranges  of  integration 
from  -1  to  1,  approximating  the  integrals  by  means  of  (18),  then  substitut¬ 
ing  (30)  for  kw  and  [1,  eqs.  (57)  and  (58)]  for  the  cosines  and  sines,  we 
obtain 


4kV*  = 


-  2j 


—  ?**>«> !  s 

K=1  L=1 


(NG) 


(NG)  2  (2)  s 

(XL  }  H1  (jINKL) 


(31) 


INKL 


where 

“INKl  *  y<xx>2  +  <yy>2  »2> 

In  (32),  XX  is  the  right-hand  side  of  (21)  with  J  replaced  by  N  and  YY  is 

the  right-hand  side  of  (22)  with  J  replaced  by  N  According  to  (28)  and 

(29),  the  XX  that  appears  explicitly  in  (31)  is  the  same  as  the  one  in  (32). 

When  J=N,  DO  loop  39  puts  the  right-hand  side  of  (31)  in  VA(I).  When  J=N, 

nested  DO  loops  40  and  41  perform  the  double  summation  in  (31) .  The  DO  loop 

indices  K  and  L  obtain,  respectively,  the  summation  indices  K  and  L  in  (31). 

DO  loop  41  accumulates  in  V2  the  summation  with  respect  to  L  in  (31).  DO 

loop  40  accumulates  in  VI  the  summation  with  respect  to  K  in  (31).  Line  80 

puts  ‘jvkl  (32)  in  XI.  Line  84  puts  the  Bessel  functions  and 

^l^INKL^  ®^^)  and  BY(2),  respectively.  Line  85  uses  (23)  with  n=l  and 

x  =  to  properly  increment  V2.  Line  95  puts  the  right-hand  side  of 

(31)  in  VA(I).  Line  105  obtains  [1,  eq .  (68)]  bv  setting  VA(N)  =  0. 

4k 

In  DO  loop  21,  —  Z  is  put  in  Z((I-1)*N+I)  where  Z ^  is  given  by 

[1,  eq.  (64)].  Thus,  as  the  index  I  traverses  its  range,  the  diagonal 

4  k 

elements  of  the  square  matrix  —  Z  of  [1,  eq.  (23)]  are  stored  in  the  singly 

4k 

dimensioned  array  Z.  Storage  of  —  Z  of  [1,  eq .  (23)1  is  bv  columns  in  the 
singly  dimensioned  array  Z.  When  I=N,  C  and  C,  of  [1,  eq .  (86)1  are  put 


in  CA  and  CP,  respectively,  and  k-'C  of  [1,  eq .  (116)]  is  put  in  CTNC. 

4k 

The  calculation  of  the  diagonal  elements  of  —  Z  of  [1,  eq .  (23)] 
is  described  as  follows.  Replacing  i  by  I  in  [1,  eq .  (64)],  multiplying 
both  sides  of  that  equation  by  4k/",  changing  the  variables  of  integration 
therein  to  make  the  ranges  of  integration  from  -1  to  1,  and  then  approxi¬ 
mating  both  integrals  by  means  of  (18) ,  we  obtain 


4k 


'II 


2 

~  •  [4  +  &  (3-2  In  (0 . 5  ' 
4 


NG 

■'  T  )  >  1  +  L 
K=1 


NG 

7  UL  ‘ 

L=1 


(33) 


where 


U1  =  a5ng) [H^2) (XX)  -  1  +  4r  ln  (0.5Y  XX)] 


(34) 


where 


XX  =  0.5', 


Y (NG)  (NG) 


(35) 


Prior  to  DO  loop  21,  line  115  puts  Y  in  GAM.  This  is  the  e'  that  appears 
in  [2,  p.  2].  Line  122  puts  in  Z1  the  square  bracketed  quantity  in  (33). 
Nested  DO  loops  22  and  23  perform  the  double  summation  in  (33).  The  DO  loop 
indices  K  and  L  obtain,  respectively,  the  summation  indices  K  and  L  in  (33). 

DO  loop  23  accumulates  in  Z2  the  summation  with  respect  to  L  in  (33).  Since 
the  square  bracketed  quantity  in  [1,  eq.  (60)]  approaches  zero  as  x  approaches 
zero,  the  term  for  which  L=K  in  (33)  is  zero  and  therefore  does  not  contri¬ 
bute  to  Z2.  Line  136  puts  XX  of  (35)  in  XX.  Line  137  puts  the  Bessel  func¬ 
tions  J  (XX)  and  Y  (XX)  in  BJ(1)  and  BY(1),  respectively.  Line  138  uses 
o  o 

(23)  with  n=0  and  x=XX  to  put  U1  of  (34)  in  U1 .  Line  139  accumulates  in  Z2 
the  summation  with  respect  to  L  in  (33).  Do  loop  22  accumulates  in  Z1 


the  quantitv  in  the 
of  (33)  in  Z( (T-l )*N+T) . 


braces  in  (33).  Line  153  puts  the  richt-hand  side 


The 

paragraph 
and  using 


calculation  of  C  and  C  s 

A 

and  the  next  paragraph. 
(18)  to  approximate  both 


of  [1,  eq .  (86)]  is  described  in  this 
Substituting  (30)  into  (1,  eq.  (89)] 
integrals  in  [1,  eq.  (89)],  we  obtain 

J  /i  -  «fc)>2  i'i 


where  U1  is  given  by  (34)  in  which  XX  is  given  V’  (35)  with  y  replaced 
by  y  •  Substitution  of  (30)  into  [1,  eq .  (96)]  gives 


IA  2  ^2  +  j(4  ln  8 


X 


))] 


( 


Substituting  (30)  and  [1,  eq.  (99)]  into  [1,  eq.  (90)]  and  using  (18)  to 
approximate  both  integrals  in  [1,  eq.  (90)],  we  obtain 


sc  A  NO  X(N« 

C,  =  j"T  +  /  K 

K= 


K 


NG  X^NG)  U1 


Vi  -  (^G))2L=i/rr^ 


(NG)  2 


where  U1  is  the  same  as  in  (36). 

When  I=N  inside  DO  loop  21,  C  of  (36)  and  C  of  (38)  are  put  in  CA 

A 

and  CP,  respectively.  Line  124  puts  I  of  (37)  in  CA.  Line  125  puts  j'1 

A 

of  (38)  in  CP.  DO  loop  23  accumulates  in  CA1  the  sum  with  respect  to  L 

in  (36).  DO  loop  23  also  accumulates  in  CPI  the  sum  with  respect  to  L  in 

(38).  Line  138  puts  U1  of  (36)  in  Ul.  Inside  DO  loop  22,  lines  149  and 

150  accumulate  C  of  (36)  and  Ct  of  (38)  in  CA  and  CP,  respec t ivelv . 

A  : 

Multiplying  [1,  eq .  (118)]  bv  k’  ,  using  (18)  to  approximate  the 
integral  therein,  and  substituting  (30)  for  kw,  we  obtain 


k-  C 


me 


NC 


me 


a;. 


(NC) 


fT- 


<x?cV 


■os (0 . 5 


v (NC)  . 
Sln 


0.5  .  cos 


When  I=X  inside  DO  loop  21,  DO  loop  22  accumulate?  in  CTN'C  the  summation 

with  respect  to  K  in  (39) .  Outside  nested  DO  Loops  21  and  22,  line  156 

_  inc  ,  .  _inc  c 

performs  the  required  multiplication  by  cos  t  so  that  k-C  ot 

(39)  will  finally  be  stored  in  CINC. 

Line  162  puts  in  CURI(I)  the  Ith  element  of  -cl100  where  iinc  is  the 

column  vector  that  satisfies  fl,  eq.  (23)].  Line  163  puts  in  CURA(I)  the  Ith 

element  of  ’'I3  where  1 3  is  the  column  vector  that  satisfies  [1,  eq .  (34)].  In 

h  q 

the  previous  two  sentences,  1=1,2,  ...,N.  Line  168  puts  2kr;Y  of  [1,  eq .  (86)] 
in  YHS .  Line  169  puts  8krYhs  in  YAB .  DO  loop  26  adds  to  YAB  the  product  of 
4k-  with  the  summation  with  respect  to  j  in  [1,  eq .  115].  Thus,  after  exit 
from  DO  loop  26,  4kD(Ya  +  Yb)  of  [1,  eq.  (115)]  will  reside  in  YAB.  Line  173 
sets  VA(N)  equal  to  -y  which,  according  to  (30),  is  2~kw  so  that  -8kCj  will 
reside  in  VA(J)  for  J=1,2,...,N.  Here  Cj  is  given  by  [1,  eqs.  (109)  and  (102)]. 
DO  loop  27  accumulates  in  TI  the  product  of  -8kr  with  the  summation  with  re¬ 
spect  to  j  in  [1,  eq .  (116)].  Line  178  puts  4k-I  of  [1,  eq .  (116)]  in  TI . 

Line  179  puts  V  of  [1,  eq .  (10)1  in  V. 

Lines  183  to  188  perform  the  additional  calculations  that  are  necessary 
to  determine  VW  of  [1,  eq.  (A-26)].  DO  loop  35  sets  VA(J)  =  0  for  J=1 , 2 , . . • ,N-1 
so  that,  after  exit  from  DO  loop  35,  ik(V^  -  )  of  [1,  eq.  (A-23) i  will  re- 

„  ,  _3W  ^bwN  rr-i 

side  in  VA(J)  for  J=1 , 2 , . . . , N .  Line  186  puts  -  ( -  1 n  ;  01  [i,  eq .  (A-24)J 
in  CURA(N)  .  Line  187  puts  4k-  (Y™  +  YbW)  of  [1,  eq.  (A-25)]  in  YABW .  Finally, 
line  188  puts  VW  of  [1,  eq .  (A-26)]  in  VW. 

Our  description  of  the  main  program  for  the  methods  of  solution  with 
and  without  pseudo-image  is  summarized  in  Table  2  where  key  variables  in  this 
program  are  listed.  Each  variable  is  identified  bv  the  line  where  it  was 
read  in,  defined  or  incremented,  and  bv  its  corre spond ing  quantity  in  [1] 


22 


Table  2. 

Key  variables 

in  the  computer  program,  program  lines  where 

they  are  read 

in,  defined  or  incremented. 

and  their  corre- 

sponding  quantities  in  [1]  or  in 

the  text 

• 

Program 

Program 

Corresponding  quantity 

Equation(s)  where 

variable 

line 

in  [1]  or  in  the 

text 

the  quantity  appears 

N 

11 

N 

[1,  eq.  (22)] 

NG 

11 

NG 

(18) 

BK 

11 

k 

[1,  eq,  (44)] 

PINC 

11 

inc 

5 

[1,  eq.  (44)] 

X(I) 

15 

XI 

[1,  eq.  (53)] 

Y(I) 

19 

yl 

[1,  eq.  (53)] 

YO 

22 

POra  and  QOm 

(7), (9), (12), (15), (14), (16) 

XG(I) 

AG  (I) 

25 

28 

(NG) 

1  (NG) 

AI 

.me 
cos  $ 

(18) 

(18) 

CSP 

32 

[1,  eq.  (44)] 

SNP 

33 

,  inc 
sin  $ 

[1,  eq.  (44)] 

DX(J) 

44 

0.5k(xJ+1  -  Xj) 

(21) 

DY(J) 

45 

°'5k(yj+i  ~ 

(22) 

G(J) 

XP(J) 

46 

47 

0  •  5'j 

[1,  eq.  (53)] 

[1,  eq.  (55)] 

YP  (J) 

48 

[1,  eq.  (56)] 

XI 

80 

'TJKL 

(20) 

BJ  (1 ) 

82 

(19),  (23) 

BY  (1 ) 

82 

Y  (arnrr) 
o  IJKL 

(19),  (23) 

Z2 

86 

The  sum  on  L 

(19) 

Z1 

88 

The  double  sum 

(19) 

Z ( (J-l ) *N+I )  93 

*  z 
n  IJ 

(19) 

Z ( (1-1 ) *N+J)  102 

*  z 

n  JI 

(24),  (19) 

SI 

96 

SI 

(26) 

VINC(I) 

100 

4kVjnr 

(25) 

VTNC(N) 

110 

4kvjnc 

(25) 

XI 

80 

aiNKL 

(32) 

BJ  (2) 

84 

Jl(  lINKL) 

(31),  (23) 

BY  (2) 

84 

Y  (<  ) 

lv  TNKL; 

(31),  (231 

Table  2.  (continued) 


Program 

variable 


Program 

line 


Corresponding  quantity  Equation(s)  where 
in  [1]  or  in  the  text  the  quantity  appears 


V2 

85 

The  sum  on  L 

(31) 

VI 

89 

The  double  sum 

(31) 

VA(I) 

95 

4kV* 

(31) 

VA(N) 

105 

4kVS 

[1,  eq.  (68)] 

GAM 

115 

Y 

(33) 

Z1 

122 

The  square  bracketed 
quantity 

(33) 

XX 

136 

XX 

(35) 

BJ(1) 

137 

J  (XX) 
o 

(34),  (23) 

BY(1) 

137 

Y  (XX) 

0 

(34),  (23) 

U1 

138 

U1 

(34) 

Z2 

139 

The  sum  on  L 

(33) 

Z1 

146 

The  quantity  in  the  {  } 
braces 

(33) 

Z((I-1)*N+I) 

153 

z 

n  II 

(33) 

CA 

124 

XA 

(37) 

CP 

125 

j71 

(38) 

CA1 

142 

The  sum  on  L 

(36) 

CPI 

143 

The  sum  on  L 

(38) 

CA 

149 

CA 

(36) 

CP 

150 

C, 

(38) 

CINC 

151 

The  sum  on  K 

(39) 

CINC 

156 

i „rinc 

knc 

(39) 

CUR I (I) 

162 

__  inc 

n  h 

[1,  eq.  (23)] 

CURA(I) 

163 

T3 

n  IT 

[1,  eq.  (34)] 

Table  2. 


(cont inued) 


Program  Program  Corresponding  quantity  Equation  (s)  where 

variable  line  in  [1]  or  in  the  text  the  quantity  appears 


YHS 

168 

2knY 

[1, 

eq. 

(86)] 

YAB 

169 

8k-Yhs 

[1, 

eq. 

(115)] 

YAB 

171 

4k- (Y3  +  Yb) 

[1, 

eq. 

(115)] 

VA(J) 

95,173 

-8kC~ 

[1, 

eqs . 

(109)  and  (102) 

TI 

176 

-8k-  multiplied  by  the 
sum  on  j 

[1, 

eq. 

(116)] 

TI 

178 

4k- 1 

[1, 

eq. 

(116)] 

V 

179 

V 

[1, 

eq. 

(10)] 

VA(J) 

173,184 

4k  (Vf  -  V$") 

[1, 

eq. 

(A— 23 ) ] 

CURA(N) 

186 

,Taw  Tbw. 

u. 

eq. 

(A-20) ] 

YABW 

187 

4k-|(YaW  +  YbW) 

[i, 

eq. 

(A-25) ] 

VW 

188 

vw 

[i. 

eq. 

(A-26) ] 

CO  1  c 
CO  2  C 
CO-3 
004 
COS 
006 
C07 
008 
009 
CIO 
01 1 
012  11 
C 1  3 
Cl  4  12 
C 1  5 
01613 
017 
018  14 
019 
02  C 
02115 
C22 
023 
02416 
02  5 
02  6 
027  17 
C2  8 
029 
030  18 
031 
C32 
033 
034 
C3  5 
C36 
037 
038 
C39 
C4Q 
041 
C42 
C43 
044 
04  5 
C46 
04  7 
04  8 
049  37 
050 
051 
C  52 
053 


LISTING  OF  THE  MAIN  PROGRAM  PCR  THE  METHODS  0? 
SOLD  IIOl!  WITH  AND  WITHOUT  PSEUDO-IMAGE 
CCMPIEX  U2,Z1,V1,Z2,72,Z(16CC),VA{40)  ,7  me  (4  0)  ,  CA 
COMPLEX  CP,CA1,CP1,U1,CURI  (4  0)  ,CUHA{40)  ,YHS,YA3,T: 
COMPLEX  V,TAEW,7W 

COMMON  10(33)  ,  PI2,  PI4,  PI7.  3  J  (1  00)  ,EY(10Q) 

DIMENSION  X  (4  0)  ,Y(40)  , XG{10)  ,AG(10)  ,  EX(4  0) 
DIMENSION  DY(40),G(4Q)  ,  XP  (  4  0)  ,?P{4  0)  ,IPS{4  0) 

OPEN  (UNIT=2  0,PI1E=’ MAUTZ1. DAT* ) 

OPEN (UNIT=21,FILE=' MAUTZ2.DAT*) 

BEAD (20, 11)  N,NG,EK,PINC 
FCSM AT  ( 213,  2E14.  7) 

WRITE  (2  1,  12)  N,NG, EE,?INC 

FORMAT ( *  N  NG* , 6X  ,  ' 3K  * , 11X, 1 PINC’/IX, 213, 2514. 7) 
READ (20,  13)  {X  (I)  ,1=1, N) 

FORMAT  (5E14.7) 

WRITE  (2  1 ,  14)  (X(I)  ,1  =  1, N) 

FORM  AT  (*  X*/ (1X,5E14.7)  ) 

READ(20, 13)  (T  (I)  ,1  =  1,  N) 

WRITE  (21,15)  (T  (I)  ,1  =  1 ,  N) 

FORM  AT  ( '  Y' /(II, -214.  7)) 

READ  (20,13)  (YO  (I)  ,1=1,33) 

WHIT E  ( 2  1 ,  16)  (YC  (I)  ,1  =  1,33) 

FORMAT  (*  YO  »/ (1X.5E14.7)  ) 

R  EAD  (20,  13)  (XG  (I)  ,  1  =  1  ,  NG) 

WRITE  (21,17)  (X G  (I)  ,1=1,  NG) 

FORM  AT ( *  XG •/ ( 1 X, 5E 1 4. 7) ) 

HEAD  (20,13)  (AG  (I)  ,1=1,  NG) 

WRIT  E  (2  1,  18)  (AG  (I)  ,  1  =  1  ,  NG) 

FORMAT  {'  AG'/ (1X,5E14.7)  ) 

EK5=. 5*EK 
CSP  =  CCS  (PINC) 

SNP=SIN (PINC) 

PI=3. 141593 
PI  2=  2. /PI 
PI4=PI/4. 

PI7=.75*PI 

02=  (0. ,-2 . ) 

DO  37  J=  1,N 
J1=J*1 

IF(J.EQ.N)  J1=1 
D  1=  8K5*  (X  (J1)  -X  (J)  ) 

D2=3K5*  (Y(JI)-Y(J)  ) 

DX  (J)  =D  1 
CY (J) =D  2 

G  (J)=SCET  (D1*D1  +  D2*C2) 

XP(J)  =BK5*{X(J1)  >X  (J)  ) 

YP (J) =  BK5* (Y  (J1) ♦?  (J) ) 

CONTINUE 
N  M=  N  -  1 
CO  38  J- 2 , N 
11=  (J-1)  *  N 
C XJ=D  X  (  J) 


C54 
055 
05  o 
C5  7 
C58 
059 
C60 
061 
062 
063 
C64 
065 
C66 
06  7 
063 
06  9 
07C 
071 
C72 
073 
C74 
075 
C76 
077 
07  8 
079 
080 
081 
082 
083 
08419 
085 
08620 
08741 
C88 
089 
C9C40 
091 
C92 
093 
C94 
09  5 
096 
097 
C9  8 
099 
100 
10  142 
102 
10339 
10438 

105 

106 


E i J=  C i  (J) 

GJ  =  G  (J) 

XPJ=XP  (J) 

YPJ=YP  (J) 

JH=J-1 

CO  39  1=1, JN 
EXI  =  DX  m 
CYI=D Y  (1) 

GI=G  (I) 

XPI=XP  (I) 

YPI=  YP  (I) 

XIJ=XPI-X?J 
YIJ= YPI-YPJ 
Z1=(0.,0.) 

IF(J.EC-N)  V1=(0.,0.) 

CO  40  K= 1 , NG 
A  K= AG  (R) 

XK=XIJ*XG  (K)  *DXI 
YK= YIJ*XG (K)  *DYI 
Z2= (C., 0.) 

IF(J.EC-N)  V  2={0.,0.) 

DO  41  L=1,NG 
XL=XG  (1) 

AL  =  AG  (L) 

X  X=  X  K-XI*DX J 
YY=YK-XL*DYJ 
X1=SCHT (XX*XX*YY*YY) 

IF(J.EQ.N)  GO  TO  19 
CALL  EES  (0,X1) 

GO  TO  20 
CALL  EES  (1 , XI) 

72=AL*SQ3T (  1.-XL*XL)/X1*CMPLX  {BJ  (2)  ,-BY(2)  )  +  7 

22=AL*C  fiPLX  (BJ  (1) ,-EY (1) ) +  Z2 

CCNTINOE 

Z1=AF*Z2*Z1 

IF(J.EQ.N)  71=AK*XX*72+71 

CONTIHOE 

GG=GI*GJ 

II«XZ*1 

Z  (II)  =GG*Z  1 

IF  (J.  IE.  8)  GO  TO  42 

VA  ( I)  =GG*U2*71 

S1=DXI*CSP*DYI*SNP 

S3=  1. 

IF  (SI  •  8E.  0.  )  S3=SIN  (S1)/S1 
S2=XPI*CSP*YPI*SNP 

7INC  (I)=8.*GI*S3*CHELX  (CCS  (S2)  ,SIN  (S  2)  ) 
JJ-(I-I)  *N+ J 
Z  (JJ)  =  Z  (II) 

CONTINUE 
CONTINUE 
VA  (N)  =0. 

GN=G  (N) 


107 
10  i 

10  . 

110 
11 1 

11229 

113 

114  30 

115 

116 

117 

118 
119 
12  C 
121 
122 

123 

124 

125 

126 

127  24 

128 

129 

130 

131 

132  25 

133 

134 

135 

136 

137 

138 

139 

140 

141 
1  42 
143 
14423 

145 

146 

147 

148 

149 

150 

151 
15222 

153 

154 
15521 

156 

157 
15828 
159 


S1=GN*SN? 

£3=1. 

IF  {Si. ME. 0.)  S3=SIN  (SI) /SI 

VINC  (N)  =  8.*GN*S3 

WHITE  (21, 30)  (VINC(I)  ,1=1, N) 

FCB8  AT ( *  VA •/ (IX, 4214. 7) ) 

WHITE  (21,  29)  (7A  (I)  ,1  =  1  ,  N) 

F0B8AT  ('  VINC */  (1X,4E14.7)  ) 

GAM  =  1.78  1072 
GA82=. 5*G AN 
IZ=  1 
NP=N*1 
CO  21  1=1,  N 
GI=G  (I) 

GA8G=GA8*GI 

Z1=C8F1X  (4.  ,4. /PI*  (3.-2.  *ALCG (GANG)  )  ) 

IF  (I.NE.N)  GO  TO  24 

CA=. 5*PI*C8  PLX ( . 5*P  I, . 25- ALOG (.25*GA8G)  ) 

CP=CBPL  I  (0.  ,PI) 

CINC=0. 

DO  22  K= 1 , NG 
Z 2=  (C.,  C.) 

IF (I. BE. N)  GO  TO  25 
CA  1=  (0.,0.) 

CP1=  (0.  ,0.) 

XGK  =  XG  (K) 

DC  23  L= 1  ,NG 
IF(L.EQ.K)  GO  TO  23 
XGL=XG  (L) 

XX=GI*ABS (XGK-XGL) 

CALL  EES  (0,XX) 

01=AG(L)  *C8PLX(BJ{1)  -1.,-BY  (1)  +PI2*AIC9  (GA32*XX)  ) 
22=01*22 

IF(I.NE.N)  GO  TO  23 
SQ=SGBT  ( 1 .- XGL*XGI) 

CA1  =  SQ*0  UCA  1 
C  P 1  =  XGL/S  Q*  U 1 +CP 1 
CCNTINOE 
AGK=  AG  (K) 

Z  1=AGK*Z2*Z  1 

IF  (I.  HE.  N)  GO  TO  22 

SQ=SGBT(1.-XGK*XGK) 

CA=  AGK*SQ*C  A 1+C  A 

CP  =  AGK*XGK/SQ*CP  1+CP 

Cl NC=CI HC  4AGK*S  G*COS (S1*XGK) 

CCNTINOE 
Z  (IZ)  =GI*GI*Z  1 
IZ=IZ+NP 
CONTINUE 

CINC=GN*CSP*CINC 
SHITE  (21,28)  (Z  (I)  ,1=1,  N) 

FCBHATC  Z*/(1X,4E14.7)) 

SHITE  (21,34)  CA,CP,CINC 


XVAIW.V.  V  ' 


"y  ■tt:  t^wv  t/ 


160  34 

161 
162 

163 

164 
16531 
166 
16732 
166 

169 

170 

171 

172  26 

173 

174 

175 

176 

177  27 

178 

179 

180 
18  1  33 
182 

183 

184 

185  35 

186 

187 

188 

189 

190  36 
9  1 
92 


F  C  3  M  A  T  { *  CA=*,2S14.7,*,  CP  =  • , 2S 14. 7/ '  CINC 
CALL  CECOMP  (N,IPS,Z) 

CALL  SCLVE(N,IPS,Z,  VTNC,C0HI) 

CALL  SOLVE (N,I?S,Z, V  A  *  COB  A) 

WRITE{21,  31)  (CUai(I)  ,1=1, N) 

PCBMAT{*  CTJSIV  (1X,4S14.7)  ) 

WRITE  (21, 32)  (CURA  (I)  ,1  =  1,  N) 

PCHMAT  (*  CORA*  /  (1X,4E14.7)  ) 

YHS=GN*GN*CA-CP 
YAB=4  . *  IRS 
CO  26  J=1,NM 
YAE=YAB*COBA  (J)*VA  {J) 

CONTINOE 
VA  (N)=2.*FI*GN 
TI= (0- , 0.) 

DO  27  J=1,N 
TI=TI*C0RI(J)  *  VA  (J) 

CONTINUE 

TI-4.*CINC-.5*TI 
V=£I/  YAB 

WRITE (21, 33)  YHS, YAP , TI, V 

FORMAT  (•  YHS=» .2E14.7, • ,  YAB=*  ,  2E 1 4. 7/ 

•  TI=*  ,2E14.7,  •  ,  V=*,2E14.7) 

DO  35  J=1,NH 
VA  (J)=0. 

CONTINUE 

CAI1  SCLVE{N,IPS,Z,VA,COHA) 

YAEW=.  5*YAB*PI*GN*CURA(N) 

VW=TI/YAB8 

WRIT E  (2 1 , 36)  TABS, VI 

FORMAT  ( *  YAEW=*  ,2E14.7,  «,  VW=*,2E14.7) 


LISTING  CF  THE  iNFOT  DATA  FILE  HAUTZ1.DAT 


36  10  0.15 
0.  C C CO 000 E 
0.  4226 1  83 E 
0. 1255014E 
0.  1902502E 
0.  19025C2E 
0.  1255014E 
0.  4226  183E 
0. OCCOOCO  E 
0.87155 74 E 
0.  819  1  520  E 
0. 9659258E 
0. 4226183E 
-0. 4226  183E 
*0.  96592  58  E 
-0. 8191 520E 
-0.  8715  574 E 
-0. 3072S82S 
0.  1600  171  E 
0.  ICCOOCOE 
-0.  1 275602E 
0.  lOOOOCOE 
0.  1829893  E 
0. 4687498E 
■0 . 97  39  065  E 
0.  1488743E 
0. 6&67134E 
0. 2955242E 


70796E 
>00  0. 
>00  0. 
♦  01  0. 
>01  0. 
>01  0. 
>01  0. 
>00  0. 
>00 
-01  0. 
>00  0. 
♦  00  0. 
>00  0. 
♦  00-0. 
>00-0. 
>00-0. 
-01 

>04  0. 
>02-0. 
♦01-0. 
>01  0. 
♦  01  0. 
-02-0. 
-01-0. 
♦00-0. 
♦  00  0. 
-01  0. 
♦  00  0. 


♦  01  0.31415! 
3026887 E- 01 
5735764E+00 
1418813E+01 
1962  121E+0  1 
1 8 1 5347E+0 1 
1Q83350E+0  1 
2  890879E+00 

2588190E+00 
9063078E ♦OO 
9063  C78E  +  00 
258  8 1 9  OE  +  O  0 
5735764E+00- 
996  1947E  +  QC- 
707 1068E+0  0* 

7366758B+04- 

596  1089E-0  5 

6024727E+04 

5832787  E- 0  4- 

9999999E+00- 

3191328E-04- 

1998720E-0  3 

865C634E+00- 

4333954E+00 

1494513E+00 

2692667E+00 


93E  +  C1 

0.8988691Z-01 
0. 7  37  37  57E  +  00 
0. 1569771E+01 
0. 1992389E+01 
0.  170  3301E+Q1 
0 . 909C  390E  +  00 
a.  1770  427E  +  00 

0.4226 1 83E+00 
0. 96  5  9  25  8E  +  CC 
0. 3  19 1  520E+00 
0.871 5574S-01 - 
-0.7071 068E+00- 
-0.9961  947E  +  00- 
-0.5735764E+00- 

-0.  6085 1 00 E  +  03 
0. 9545773 E- 08 
0.  1613512E  +  04- 
-0. 1107698E-06 
-0.  1097659E-02 
-0. 1562498E-01 
0.7317495E-05 
-  0. 6794  096E  +  00- 
0 .6794C96E  *00 
0.2 190864E+00 
0. 2190864E+00 


0. 1770427E+00  0. 
C.9090390E+00  0. 
0.  17033C1E+Q1  0. 
C. 1 992 389E*  0 1  0. 
0. 1569771E+0 1  0. 
C.7373757E+00  0. 
0. 898  86  9 1 E-0  1  0, 

0. 57357 64E+00  0. 
C.9961 947E+00  0. 
0. 7071 068E+90  0. 
■C.8715574E-01-0. 
■C. 81915  20E+ 0  0-0. 
■0. 9659258  E+ 00-0. 
■0.4226183E+00-0. 

0.  1710234E  +  02-0. 
0.4 1631 50E+05  0. 
■0.  753  22  10E  +  Q  2  0. 
G.3072946E+05  0. 
0.2461 455E-04  0. 
0.1427079E-03-0. 


28908  79E+  00 
,  1083350E  +  01 
,  18  15347Z01 
,  196212  1 E  +  0 1 
,  14  188 1 3E  +  01 
,  5735764  E  +  00 
,  3026887E- 01 

,  7071068E+00 
,  996 1 947E  +  00 
,  5735764E+00 
2588 190E+00 
,  9063 078E+  OO 
906  3078E  +  00 
258819CE+00 

,  227100  1E  +  00 
,  342021  1E+03 
.  1402590E  +  01 
,  288643  1E+03 
,  1000000E+C1 
59374  34E-05 


0.4333954E+00-0. 1488743E+00 
0. 8650634E+00  0. 973906 53+00 
0.2692667E+00  0. 2955242S+00 
0.  1494513E+00  0 . 6667 1 34E- 0 1 
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C 

C 

N 

36 

X 

0. 

0. 

0. 

0. 

0. 

0. 

0. 

0. 

Y 

0. 

0. 

0. 

0. 

-0. 

-0. 

“0. 

-0. 

TO 

-0. 

0. 

0. 

-0. 

0. 

0. 

0. 

XG 


LISTING  CF  THZ  CCTPCT  DATA  EIIZ  MAUTZ2.DAT 

NG  BK  EINC 

10  0.  1570796E  +  0  1  0. 3 1 41 593E+C 1 


COOOOCOE+OO 
42261832+00 
12550142+01 
1902502E  +  0  1 
19025022  +  01 
12550  T42  +  01 
42  26 1 83  E  +  00 
CCCCCCOE+OO 


0. 3026887 E- 01 
0-57357642+00 
0- 14188132+01 
0. 1962 121E+0 1 
0.  18 1 5347E  +  0 1 
0. 1083350E+01 
0-289C879E+00 


87155742-01  0. 
819  1  520  E+  00  0. 
96592582+00  0. 
42261832+00  0. 
42261 83  E  +  0  0-0. 
96  5  9258E  +  00-  0. 
8191  5202  +  00-0. 
87 15574  E-0  1 

307  2582  E  +  0  4  0. 
1600171  E-0  2-0. 
1000000E  +  0  1-0, 
12756C2E-01  0. 
1000000  E  +  01  0. 
1 8 29 8 93 E— 02-0. 
4687498 E- 01-0. 


2588 1 9  OE+O  0 
9063  078E  +  00 
9063078E+00 
258ei90E+00 
57357642+00* 
996  1947E+00- 
707 1 C68E+00- 


0. 8988691E-01 
0. 7373757E+00 
0.1569771 E+01 
0.  1992  389E  +  0 1 
0-  1703301E+01 
0.9090390E+00 
0.  1770427E  +  00 


0.4226 183E+00 
0.9659258E+00 
0.8191520E+00 
0.8715574  E— 01- 
0.70710682  +  00- 
0. 996 1 947E+0Q 
0.5735764E+00- 


7368758E  +  0  4- 
59610892-05 
6024727E+04 
58327972-04- 
9S99999E+00' 
31913282-04- 
1S9872CE-03 


■0.6085  1 0OE  +  O3 
0.9545773E-08 
0.  161351 2E+  04* 
■0. 1107698E-06 
•0.  1C97659E-02 
■0.  1562498E-01 
0.7317495E-05 


0.  1770427E  +  00  0. 
0.9090390E  +  00  0, 
0.  17033C1E+01  0, 
0. 1992389E+01  0 
C.1569771E+01  0 
0.7373757E+00  0 
C. 898  86 91 E-0  1  0 


0.  5735764E  +  00  0 
C. 996 1 9  47E+  0  0  0 
C.7071G68E+00  0 
-C.8715574E-01-0 
-C.3191520E+00-0 
-C.9659258E+00-0 
■0.42261 83E+00-0 


0.171 0234 E+ 02-0 
0.41631  50 E  +  05  0 
0.75322  10E+02  0 
0.  30729  46E+05  0. 
0.2461 455E-04  0. 
0.  1 4 270  79E- 0  3-0 


2890879  E  +  00 
,  1083  35CE+01 
1815347E+01 
19621 2  1E  +  01 
141881 3E+0 1 
.  5735764E  +  00 
.  3026887E-0 1 


7071068E+00 
,  996 1 947E  +  00 
,  5735764E  +  00 
2583 190E+00 
.  9063078E  +  C0 
.9063078E+00 
,  258819CE  +  00 


,  227100 1E+00 
,  34202  1  1E  +  03 
1402590 E+01 
2886  4  3 1 E  +  33 
1 000000  E  +  0 1 
59  374  3  4E-05 


-0. 9739065E+0 0-0. 86 5C634E+00—0.6794096E+0 0-0. 43339 54 E+00 -0. 148874 3 *+00 
0. 1488743E+00  0.4333954E+00  0. 6794096E+00  0. 86 506 34 E+ 00  0 . 97 390 6 5E+ 0 0 
AG 

0. 6667 1 34E-0 1  0. 1494513E+00  0. 2 190864E +00  0. 2692667E+00  0 . 29552422+00 
0. 2955242E+00  0. 2692667E+00  0.2 19 0864E+00  0. 1 49 45 1 3E+ 00  0. 6667 1 34 E-0 1 
VINC 

0. 1094818E+01-0.26C3224E-01  0. 1089 96 OE +01-0. 1031663 E+00 
0.  1070  4  14E  +  01-0.2277548E  +  00  0 .  102 1 332 E  +  0 1- 0. 39 1 5 596E  +  0 0 
0.9268360E+00-0.5797543E+00  0.7749375E+00-C.7703197E+00 
0.  5626 479E  + 00- 0.936 139  IE +0  0  0 . 2993  4  8  0E +00-0.  10 500 82E+ 0  1 
0. 6 526236 E- 02-0. 109 179 1E+0 1-0. 2867749E+00- 0. 1 053589E+ 0 1 
-0. 55 14 164E+ 00- 0. 942 7966E+00-0. 7656 7 32E+00-C. 7795281 E+00 
-0. 9 19 8 373 E+00- 0. 5 9079 19E +00-0. 1 C 16 58 1 E+ 0 1- 0. 4037424E+ 00 
-0.  1067612E  +  01-0.2405345E+00-0.  1038 65 1 E  +  0  1-0.  1161893E  +  00 
-0.  109 4 428 E  +  0 1-0.  39 1  1 87 2 E- 0  1- 0 . 1 09 5 1 5 3E  +  0 1 - C.  1 3 09 4 29 E-0  1 
-0.  109  4  42  8E+01-0. 391  1890E-01-0.  1088651 5  +  0 1-0.  1  16  1897E+ 10 
-0.  1067 6 12E +  01-0. 240 5 35  IE +  00-0.  10 16  5  8 1E  + 01- C. 40 374  31 E+00 
-0.9198363E+00-0.5907926E+00-0.7656725E+00-0.7795288E+CO 
-0.  55 141 55 E  +  00- C. 9427971 E+0 0-0. 28 6  77 39E  + 00- 0.1  0535 89E+ 01 
0.  6527  342E-02-0.  109 1791 E  +  01  0. 299349 1E  +  00-0.  10  5008  1E+ 3  1 
0. 5626488E+00-0.9361386E+00  0.7749382E+00-C.7703190E+00 
0.  9268364S  +  00-0.579753  5E  +  0  0  0 .  102  1  332E  +  0  1-0.  39  15 5 90E+ 3 0 


0.  10704l4E*01-0.2277593E*00  C.  1C8996 0E  +  01-0.  1031659 E*  00 
0.  1094818E*01-0.26032C4E-0  1  0.  1 0 9 5 23 1 E +  0  1  C . C 0 00 C OOE* 0 0 
7  A 

0.2994785E-01-0.1381276E-02  0. 2858944 3-0 1-0. 5333051 E- 02 
0.  2962683E-0 1-0.  1  1316  ICE- 0  1  0 . 2 9 1 1  25 8E-0 1 - G. 1 85456 1 E-0  1 
0. 2649 176E- 01 -0.26 13214E-01  0. 2 1 76669E-0 1- 0. 33229 22E-0 1 
0.  15281 27 E- 01-0. 3915687 E- 01  0. 7584020E-02-0.4347779E-0  1 
-  0. 69  80846E-03-0 . 4602249E-0 1-0.8966  518E-02-C.4  6866  16 E-0  1 
-0.  1672587 E—  0 1-0.4627063  E-0 1-0. 2362260 E- 01-0. 4461700E-0 1 
-0. 29449 13E-0 1-0. 4233324E-0 1-0. 34 12371E-01-C.3983858E-01 
-0. 37 65389E-01-0. 37 50189 E- 01-0. 4C09859E-01-0. 3562038 E- 01 
-0.  4152818E-0  1—0. 3440646 E-0  1- 0. 4199 799E-01-C. 33987 86 E-0 1 
-0.  4152818E-01-0.3440646E-01-0.4009  859E-01-0.3562038E-0  1 
-0. 3765389E-01-0.375C189E-01-0.3412371E-01-C.3983858E-01 
-0. 29449 13E- 01- 0.4233324 E-0 1-0.  23622 60 E-0  1-0. 4 46  17 00 E-0  1 
-0.  1672587E-01-0.4627G63E- 01-0. 8966518E-02-C. 46866  16E-01 
-0.6980  846E-03-0.4602249E-01  0. 758  4  020E-0  2-0.  4  347779E-0  1 
0.  15281 27 E- 01-0. 3915687 E- 01  0.2 176669E-01-0.3322932E-01 
0.  264  5 176E- 01-0.26 13 2 14 E-0  1  0. 2911 258E-0 1-0.  185456 1E-0  1 
0.  2962683 E-0 1-0.  113 16  10 E-0  1  0. 2 85 8944E-01 -0. 53 33 05 1 E-0 2 
0.2594785E-01-0- 138 1276 E- 02  0 . OOOOOOOE *00  C. COOOOOOE+O 0 
Z 

0.  7473691E-01  0. 13850262*00  0.7335540E-0 1  0.7065693E-0  1 
0. 6932754 E-01  0. 2968639E-0 1  0. 6298 721E- 0 1  C.6236447E-02 
0.  5483766E-0 1-0. 1027953 E-01  0. 45484 06E-0 1-0. 22 16  265E-Q  1 
0.  3 5 55794 E-0  1-0. 303 74  57 E-01  0 . 2564739E-01 -0. 3 554 059 E-0 1 
0.  1624 3 63 E-0 1-0. 382 0945 E-01  0. 7709407E-02-0. 38 9 1708 E-0  1 
0.  2694  7 17  E- 03- 0. 38 18 3 08 E-0 1-0. 597 8546 E-0 2- 0.36 492 54E-01 
-0.  1  104134E-01-0.3427249E-0 1-0. 1499958E-01-0. 3187981  E-0  1 
-0.  17  97623 E- 0 1-0- 2959  58  1E-0 1-0. 2C 1C  92 6E- 01- C.276 3 171 E-01 
-0.2152729E-01-0.2613496E-01-0.2233387E-01-0.2520 174E-0 1 
-0. 2259513E-01-0.2486496E-01-0.2233394E-C1-C.2520178E-01 
-0.  21 52726E- 01- 0.26 13493E-0 1-0. 2  010 927E-0 1-0. 27631 73E-0 1 
-0.  1797624 E-0 1-0. 295 9 5 84 E-0 1-0. 1499553E-0  1-C. 31 87969E-0  1 
-0. 1104136E-01-0.3427249E-01-0.5978571E-02-0.3649271E-0 1 
0. 2694886 E-03- 0.38  18 30 3E- 01  0. 770 9 373 E- 02- C. 3 8 5 1 6 93 E-0 1 
0.  16 24 36 7E -01-0. 3820957 E-01  0. 2564739E-0 1-0. 3554058E-0  1 
0. 35557S3E-0 1-0. 30 37456 E-01  0. 454 8407E-0 1 -  0. 22 1 6266E-0 1 
0. 54837 €5 E- 01-0. 1027953 E-01  0. 6298724E-0  1  0. 6236455E-02 
0. 6932755E-0 1  0.2S68641E-01  0. 7335 536E-01  C. 70656 92E-0  1 
CA=  0.246161 6E*01  0 . 4775566 E*0 1 ,  C?=  0.  184 1 846E-0  1  0. 3 1 71688E*0 1 
CINC=-0.2151464E*00 
COBI 

0.2 11 3638 E* 01-0. 536 C954E*00  0 . 1 99 8972E+01- 0. 6872047E* 0 0 
0.  1 79 1 4  95E* 0 1-0. 903844  5 E  +  0  0  0.1479588E*01-0.  1134857E*01 
0.  10685 54 E* 0  1-0.  13 16  47 6E*01  0 . 59 17  22 8E* 00- C.  1 3869 05 E* 0  1 
0.11 102C2E+00-0.  1306445E+0  1-0. 2 97 1  143E*00-0.  1075  1423*01 
-0.  5667900 E *00-0. 7376  5372*0  0-0.  66  7178  IE* 0  0- 0. 36  97278 E* 0 0 
-0.  6132755E*CG-0.5098838E-01-0.4569358E*00  0.  16287 54 E* 0  0 
-0. 2636858E*00  0 . 25 5 20 1 4 E+00- C . 8820872E-0 1  0.24575962*00 
0.3S96931E-01  0.  1 756236E*00  0.  1 1 67466E*00  C. 8950505E-0  1 
0.15346652*00  0. 23 29606 E-0 1  0. 1 637747E+C0-0. 1 162453E-02 
0.  15346652*00  0.2329659E-0  1  0.  1 16 7 467E *0 0  C. 8950505E-0  1 


0. 3996842E-01  0 .  175623 9 S  +  00- C . 38 209 1  1E-0 1  0.  2457597E  +  00 
-0. 2636863E+00  0. 255 20 1 2E+00-0. 4 56 9 371 E *00  C.  1628750E+03 
-0. 6122758E*C0-C.50S8913E-01-0.667173lE+CQ-C.2697297E+00 
-0. 5667889E+00-0.7376546E+00-0.297 1131 E +00-0. 1075 144E+0 1 
0.  1 1 102  17E  +  00-0. 13G6445E«-01  0. 5  91 72  45E+G0- C.  1  3869  05E  +  O  1 
0.  1C6  85  55E  +  01-0.  13 1 6475E+0 1  0. 14  79589E  +  0 1-0.  1 134857E+0 1 
0.  1791496E*01-0.9038439E+00  0.  1998573E+01-G.6872036E+00 
0.  2113638E  +  01-0.536C95mO0  0. 2  14982  4E  +  0 1-C.  48  220  35E  +  0  0 
COB  A 

0. 8435730E-01-0.2000979E+00  0. 76896 5 9E-0 1-0. 74 14376E-02 
0.6570814E— 01-Q.487CQ10E-01  0 . 52363 92E-01 -  0. 51 36 077 E-0 1 
0.  3856174E-0 1-0. 553481  IE-01  0.2580129E-01-C.5480544E-0  1 
0.  1514947E-01-G.5 149234E-01  0. 7 139 C7 0E- 02- C. 4639276E-0 1 
0.  1799082E- 02-0.4062632E-0 1-0. 12190 58E-02-0.  35076  16E-01 
-0. 24729C5E-02-0.3032 159E-0 1-0. 25681 23E-02- 0. 2663 050E-0 1 
-0. 20 3 5613 E- 02-0. 240  1857E-01-0.  128 8466E-02- 0.22336 53E-01 
-0. 56944 95 E-03-0.2135902 E-0 1-0. 1 499 7 1 9E-04-0. 208 54 4 1E-0 1 
0.  3240  8C4E-03- 0. 206  3298E-0  1  0.4  37 1 C 95E-03- C.2 0573 2 1 E-0 1 
0. 324C795E-03- 0.2C63297 E-0 1-0. 1499475E-04-0. 2085441 E- 01 
-0.  56945 95 E-03-0. 2135902 E- 01-0.  1 2 6 8 466E- 02- 0.22336 53E-0 1 
-0. 203 96  18E -02- 0-2401 856 E-01-0. 2568 109E -02-0. 2663050 E- 31 
-0. 2472912E-02-0.3C32158E-01-0. 1 2 1 9C56E-02-0.35 076 17E-01 
0.  1799C82E-02-0.4062631E-01  0. 71 39078E-02-0.  4639276E-0 1 
0.  15 149 47 E- 01-0. 51 49233 E- 01  0. 258 0 1 27E-01-0. 5480543E-0 1 
0. 38  56 175E-0 1-0. 553  481  IE-01  0. 5236  392E-0 1-0. 5 136C79E-0 1 
0.6570815E-01-0.4870008E-01  0.7689659S-01-0.7414370E-02 
0. 8435730E-01-0.2CGC979E+00  0. 869780 1E-0 1  0. 3527458E+00 
YHS-=  0.  277188GE-G1-0.2082181E+C1,  YAB  =  0. 9 1  96 463E-0 1 -0 . 1 23  52  1 8E  +  02 
TI  =  -0.  1835366E+01  0. 4  1554 12E+00 ,  V=-0.3474545E-0 1-0.  1493277E  +  00 
YABi=  0. 8403454 E-0 1-0. 107705 1E+ 02,  VU=- 0. 3 9 90850E-01 -0 . 1 70095 1 E+00 
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V.  THE  SUBROUTINE  BESJ1 

The  subroutine  BESJ1(N,  X,  BJ1)  puts  J  (I*X)  in  BJ1(I)  for  1=1,2,...,: 
Here,  is  the  Bessel  function  of  the  first  kind  of  order  one,  N  is  an  in¬ 
teger  not  less  than  one,  and  X  is  a  non-negative  real  number.  The  minimum 
allocation  for  BJ1  is  given  by 

DIMENSION  BJ1(N) 

The  subroutine  BESJ1  is  listed  at  the  end  of  this  section.  If  I*X£3, 

then  BJ1(I)  is  calculated  by  lines  7  to  20  inside  DO  loop  15  of  the  listing 

As  suggested  in  [2,  Sec.  9.12.,  Example  1],  JMZ(XX)  and  (XX)  are  set 

-20 

equal  to  the  arbitrary  values  of  zero  and  10  ,  respectively,  where 

XX  =  I*X  (40) 

and 

MZ  =  11  +  2[XX]  (41) 

where  [XXl  is  the  largest  integer  that  does  not  exceed  XX.  This  is  done 

in  lines  7  to  9.  In  DO  loop  16,  the  recurrence  relation  (4)  is  used  to 

calculate  '  J  (XX),  n  =  MZ-2,  MZ-3, - ,0;.  The  calculated  value  of  J,  (XX) 

n  i 

has  to  be  divided  by  ot  of  (5).  This  division  is  done  in  line  20. 

If  I*X>3,  then  BJ1(I)  is  calculated  by  lines  22  to  32  inside  DO 
loop  15.  The  polynomial  approximation  [2,  eq.  (9.4.6)]  is  used.  With 
regard  to  [2,  eq.  (9.4.6)],  lines  28  and  29  put  f  in  F,  lines  30  and  31 
put  r  in  T,  and  line  32  puts  J  (I*X)  in  BJ1(I). 


C01C  LISTING  OF  THE  SUESCUTINE  EES  J 1 

GO 2  SUE ROUTINE  E ES J 1  { N , X , 3 J 1 ) 

CG3  DIMENSION  3 J  { 1  8 )  ,  EJ  1  ( 1  0 3 C 0) 

CQ4  DO  15  1=1, N 

COS  XX=FICAT  (I)  *X 

006  IF  { XI-  3.  )  17,17,18 

CO  7 17  EZ=1 1  +  2*IFIX  (XX) 

008  B  J  ( MZ ♦  1)  =0. 

009  B J (HZ) = 1- £-20 

010  H  1=MZ-  1 

Oil  X2=2./XX 

012  DO  16  K=1,81 

Cl  3  MK=  flZ-K 

014  EJ  (M K)  =X2*FLOAT  (flK)  *E J ( MK* 1 ) -E J  (  HK  +  2) 

01516  CONTINUE 

016  ALP=.  5*B  J  ( 1) 

017  DO  19  J=3  , HZ , 2 

018  ALP  =  AIP  *BJ ( J) 

C 1 9 19  CONTINUE 

020  8  J  1  (I)=BJ  (2)  /  (2.*ALP) 

021  GO  TC  15 

02218  X 1=  3./XX 

C23  X2=X1*X1 

024  X3=  X2*X  1 

C25  X4—  X3*  X 1 

026  X5=X4*X 1 

027  X6-=  X5*X  1 

02  8  F  =  .  7978  8 46*.  156E-05*X1* .  1 65 9667E-0 1* X2 *.  17  1  05 E-0  3*X 3- 

C29  1  .249511E-02*X4+. 1 1 3653E-02*X5-. 200 33E-0 3*X6 

030  T=XX-2.  3  561  94+.  124996  1*X1«-.565E-04*X2-.  637879 S-02*X3* 

C31  1  .74348 E-03*I4+.79824Z-03*X5-. 29166 E-Q 3* X6 

032  EJ  1  (I)  =F/SQ3T  (XX)  *CCS  (T) 

03315  CONTINUE 

C  34  RETURN 


VI.  THE  SUBROUTINE  BESJY 


The  subroutine  BESJY (Nl,  \’2,  X,  BJY)  puts  J^.CX)  Y  (X)  in  BJY(N)  for 

X=N1,  Nl+1,  N2  where  and  Y^T  are,  respectively,  the  Bessel  functions 

of  the  first  and  second  kind  of  order  N.  Moreover,  X  is  a  non-negative  real 

number,  N 2  _>  Nl ,  and  Nl  is  an  integer  appreciably  larger  than  X  so  that  Debve's 

asymptotic  expansions  [2,  eqs.  (9.3.7)  and  (9.3.8)]  apply  to  (X)  and 

Y  (X) .  The  minimum  allocation  for  BJY  is  given  by 
N 1 

DIMENSION  BJY (N2) 

Combining  [2,  eqs.  (9.3.7)  to  (9.3.9)],  we  obtain 

4  u.(t)  4  (-l)V  (t) 

J  (X)Y  (X)  =  -  4  (1  +  l  -^— )  (1  +  l - r1 - )  (42) 

-N  ‘  i=l  N  i=l  N 


where 


and 

uL(t)  =  (3t  -  5t3)/24 

u2(t)  =  (81 t2  -  462t4  +  385t6) /1152 

u  (t)  =  (30375t3  -  369603t5  +  765765t/  -  425425t9) /414720 
u,(t)  =  (4465125t4  -  94121h7hth  +  349922430tR 

4 

ID  1  > 

-  446185740t  +  1 8 39 1 072 3 1  “) / 398 131  20 

Equation  (43)  is  verit  ied  in  tin*  following  manner.  From 
f  2 ,  eq.  (9.3.7)], 


(43) 

(44) 

(45) 

(46) 

(47) 


Subst ituting 


sech  :<  =  X/N  (49) 

into  the  identity  [2,  eq.  (4.5.17)] 

tanh  :»  =  \Jl  -  sech^  (50) 

and  noting  that  coth  a  is  the  reciprocal  of  tanh  a,  we  obtain 

,  1 

coth  a  =  ~  (51) 

y i  -  (x/n)2 

Substitution  of  (51)  into  (48)  gives  the  desired  result  (43). 

The  subroutine  BESJY  is  listed  at  the  end  of  this  section.  The  right- 

hand  side  of  (42)  is  calculated  inside  DO  loop  11  of  this  listing.  The 

order  N  that  appears  in  (42)  is  obtained  as  the  index  of  DO  loop  11.  Lines  5 
2  3  4 

to  8  put  N,  N  ,  N  ,  and  N  in  FI ,  F2,  F3,  and  F4,  respectively.  With  regard 

to  (43)-(47),  lines  10  to  20  put  t1  in  TI  for  1=1,2,3,4,5,6,7,8,9,10,12. 

2 

Line  21  puts  u^(t)/N  in  Ul,  line  22  puts  u?(t)/N  in  U2,  line  23  puts  u^(t)/N 
in  U3 ,  and  lines  24  and  25  put  u^(t)/N^  in  U4.  Line  28  puts  the  right-hand 


side  of  (42)  in  BJY(N). 
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001C 
C02 
C03 
CO  4 
00  5 
C06 
007 
COS 
C09 
010 
01 1 
Cl  2 
013 
Cl4* 
015 
Cl  6 
017 
Cl  8 
019 
020 
021 
022 
023 

024 
025 
C2  6 
027 
02  8 
029  11 
C30 
031 


LISTING  OF  THE  SOBBCUTINE  BESOT 
SOBFCOTINE  EES J Y { N 1 , N2, X, B JY) 

DIMENSION  B JY  ( 1  COCO) 

CO  11  N=N1,N2 
F  1  =  N 

F2=F 1*F  1 
F  3=F  2  *F  1 
F4=F3*F 1 
XN=X/F 1 

T1  =  1./SCBT(1.-XN*XN) 

T  2=T  1  *T  1 

T3=T2*T1 

T  4=T  3*T  1 

T5=T4*T1 

T6=T5*T  1 

T7=T6*T1 

T8=T7*T1 

19=18*T1 

T10=T9*T1 

T12=T10*T2 

0  1-=  (3.  *T  1-5.  *T3)  /(24«*F1) 

02= (81.*T2-462.*T4+38  5.*T6)/  (11S2.*F2) 

03= ( 30375. *T3-369603.*T5+765765.*T7-425425. *T9) / (414720- *F3) 
04=  (44651 25. *T4-94 1 21 68. E+0 1*T6+3499224. E+02*T8-446 1 857. E+02* 
T10+18591  07.E+  G2*T  12)/(  39  81312.  E+01*F4) 

05=1.+U2*04 

06=01+03 

BJY(N)  =  -T1/<3.  141593*F1) *(05+06) *(05-06) 

CONTINUE 

BETOHN 

END 


* 


m 


.•  k  1  .  a  _  a  _  a 


THE  MAIN  PROGRAM  FOR  THE  FOURIER  SERIES  METHOD  OF  SOLUTION 


The  main  program  for  the  Fourier  series  method  of  solution  calculates 
the  complex  constant  V  of  [1,  eq.  (B— 25) ] .  Input  data  are  read  from  the 
file  MAUTZ3.DAT,  output  data  are  written  on  the  file  MAUTZ4.DAT,  and  the 
subroutines  BES,  BESJ1,  and  BESJY  are  called.  These  subroutines  are  de¬ 
scribed  in  Sections  II,  V,  and  VI. 

The  main  program  for  the  Fourier  series  method  of  solution  is  listed 
at  the  end  of  this  section.  In  the  open  statements  on  lines  6  and  7  of  this 
listing,  MAUTZ3.DAT  is  the  input  data  file  and  MAUTZ4.DAT  is  the  output  data 
file.  The  input  data  are  read  early  in  the  program  according  to 

READ (20, 17)  Nl,  N2,  X,  P 
17  FORMAT (13,  15,  2E14.7) 

READ(20,11) (Y0(I) ,  1=1,  33) 

11  FORMAT (5E14. 7) 

The  definitions  of  Nl  and  N2  are  based  on  the  technical  digression  in  the 
next  paragraph. 

Truncating  the  infinite  series  in  [1,  eq.  (B— 2 5 ) 1 ,  we  obtain 

2  m  "jka  cos 

»  -  *  f§  *  <52) 

o 


where 


.  Nl-1  .N  J. (Nl  ) 

U1  =  - TJT -  +  l  -TtJ -  (  1  .  ~  ) 

4H^2)(ka)  N=1  H^2)(ka)  NCo 


J,  (ND  ) 


Jo(ka)H^'(ka)  N=1  .^(ka)^"'  (ka) 


N2  ,  J, (Nl  ) 

—  _ 1 _  ,_1 _ o_ 

+  J  -  .  I. (ka)Y.. (ka)  ^  NO 


is  the  Hankei  function  of  the  second  kind  of  order  N 


In  (53)  and  (54),  H,J2) 

According  to  (23), 

H!J2)(ka)  =  JN(ka)  -  jYN(ka)  (55) 

where  J  and  are,  respectively,  the  Bessel  functions  of  the  first  and 
second  kinds  of  order  N.  Still  in  (53)  and  (54),  N1  and  N2  are  positive 
integers  such  that 


N2 

2  Nl 

(56) 

Nl 

>>  ka 

(57) 

N2 

»  m/0 

0 

(58) 

The  inequality  (57)  allows  the  series  in  (53)  to  be  truncated  at  N=N1-1. 

It  allows  deletion  of  the  quantity  J^(ka)  on  the  right-hand  side  of  (55) 

(2) 

whenever  N  2.  Nl«  so  that  j/YN(ka)  instead  of  l/l^  (ka)  could  be  put  in 

the  second  summation  on  the  right-hand  side  of  (54).  The  inequality  (58) 

is  necessary  because,  in  order  for  [1,  eq.  (B— 9 ) ]  to  be  well-satisfied, 

e^m*  terms  whose  periods  are  appreciably  smaller  than  2$o  must  be  retained 

on  the  left-hand  side  of  [1,  eq.  (B— 9)  ]  .  The  period  of  e'5™'^  is  20 when 

m  is  tt/v>  . 

o 

The  Nl  and  N2  that  appear  in  the  first  read  statement  are  the  same 
as  those  in  (53)  and  (54).  In  the  same  read  statement,  X  is  the  electrical 
length  ka  in  (53)  and  (54),  and  P  is  the  angle  0q  in  (53)  and  (54).  Here, 

P  is  in  radians.  In  the  second  read  statement,  the  array  YO  contains  input 
data  for  the  subroutine  BES.  The  values  of  the  elements  of  YO  were  given 
in  Section  II.  These  values  should  not  deeply  concern  the  user  because  he 
will  never  have  to  change  them. 


Minimum  allocations  are  given  by 


COMMON  BJ(N1),  BY(N1) 

DIMENSION  BJ1 (N2) ,  BJY(N2) 

Immediately  after  the  main  program  at  the  end  of  this  section, 
the  contents  of  the  input  data  file  MAUTZ3.DAT  and  the  output  data  file 
MAUTZ4.DAT  are  listed  when 

N1  =  20  ' 

N2  =  10000  (  (59 

X  =  t/2 

P  =  tt/36 

The  values  of  X  and  P  in  (59)  obtain  the  example  of  [1,  Section  VIII]. 
These  values  mean  that 

ka  =  tt/2  (60 

and 

-7  5  -  36  (61 

o 


In  view  of  (60)  and  (61),  the  values  of  N1  and  N2  in  (59)  satisfy  the 

inequalities  (57)  and  (58).  The  output  data  file  MAUTZ4.DAT  contains  all 

the  data  put  out  by  the  write  statements  in  the  main  program  for  the 

Fourier  series  method  of  solution.  The  contents  of  the  output  data  file 

MAUTZ4.DAT  are  described  in  the  next  three  paragraphs. 

The  input  data  are  written  out  immediately  after  they  are  read  in. 

The  ith  number  in  the  ith  row  under  the  heading  "BJ"  is  Jr/.  .  .  . .  .(ka). 

J  5(i-l)+j-l 

The  jth  number  in  the  ith  row  under  the  heading  "BY"  is  Y^^^^  ^(ka). 

By  means  of  (55),  the  Bessel  functions  mentioned  in  the  previous  two 

(2)  (2) 

sentences  are  used  to  construct  the  H'  (ka)  and  the  (ka)  that  appear 


in  (53).  The  number  preceded  bv  "B.r'_(l)  =  "  is  J,  ( :  ) ,  and  the  number  pre- 

I  o  r 

ceded  bv  ”3J1(N2)="  is  J,  (N2*T  ). 

L  0 

Each  line  under  che  heading  "N  L’i  U2"  contains  five  numbers.  The  second 
and  third  numbers  are  the  real  and  imaginary  parts  of  the  right-hand  side  of 
(53)  with  hi— 1  replaced  by  the  first  number.  The  fourth  and  fifth  numbers  are 
the  real  and  imaginary  parts  of  the  right-hand  side  of  (54)  with  the  maximum 
value  of  N  therein  equal  to  the  first  number. 


The  number  preceded  by  "BJY(N1)="  is  (ka)Y^  (ka)  ,  and  the  number 

preceded  bv  "BJY(N2)="  is  J  _ (ka)Y  _ (ka) .  The  numbers  written  under  the 

N2  N2 

heading  "BJ1"  are 


I 


I  i  MN*>2 

S-Nl  Vk3)  YN(ka)  (  N*0  }  (62) 

where  I  runs  from  N1  to  N2  in  steps  of  50.  These  numbers  are  written  in 
order  to  observe  the  rate  at  which  (62)  converges  as  I  increases.  Finally, 
the  two  numbers  preceded  by  "V="  are,  respectively,  the  real  and  imaginary 
parts  of  V  of  (52)  . 

Lines  17  to  19  define  the  common  variables  (2)  that  are  input  data 

for  the  subroutine  BES.  Line  21  puts  J^(ka)  in  BJ(N+1)  for  N=0,1, . . . ,N1-1 . 

Line  21  also  puts  Y^(ka)  in  BY(N+1)  for  N=0 , 1 , . . . ,N1-1 .  Line  26  puts 

J. (NO  )  in  BJ1(N)  for  N=l,2 . N2. 

1  o 


Line  30  puts  in  U1  the  term  not  governed  by  the  summation  sign  in  (53). 
Line  31  puts  in  U2  the  term  not  governed  by  either  one  of  the  summation 
signs  in  (54).  Inside  DO  loop  13,  line  47  adds  to  U1  the  term  whose  summa¬ 
tion  index  is  N  in  (53),  and  line  48  adds  to  U2  the  terra  whose  summation 


index  is  N  in  (54) . 


J, (NO  ) 

Line  39  puts  (— — — )  in  B. 


Line  41  puts  J  (ka)  in 


M.f. 


BJN. 


Lines  42  to  45  put  ( 


VNV 


NC 


)/R^2)(ka)  in  U3 . 


If  '  J  (ka)  '  _  10 


-10 


o 

(2) 

then,  as  a  measure  to  avoid  an  underflow,  (ka)  is  approximated  by 

N 

-jY^(ka).  After  line  46  is  executed,  jl  will  reside  in  U. 

Line  51  puts  JN(ka)  YN(ka)  in  BJY(N)  for  N=N1,  Nl+1 , . . . ,  N2.  Inside 
DO  loop  14,  line  58  accumulates  in  S  the  second  summation  on  the  right- 
hand  side  of  (54) .  The  index  N  of  DO  loop  14  obtains  the  summation  index 

Jl(NV 

N  of  (54).  Line  57  puts  ( — ^ - -  )  in  B.  Line  63  puts  the  right-hand  side 

o 

of  (54)  in  U2.  Finally,  line  65  puts  V  of  (52)  in  V. 

Our  description  of  the  main  program  for  the  Fourier  series  method  of 


solution  is  summarized  in  Table  3  where  key  variables  in  this  program  are 
listed.  Each  variable  is  identified  by  the  line  where  it  was  read  in,  de¬ 
fined,  or  incremented,  and  by  its  corresponding  quantity  in  the  text. 


Program 

variables 


Program 

Program 

Corresponding  quantity 

Equation (s)  where  the 

variables 

line 

in  the  text 

quantity  appears 

Nl 

8 

Nl 

(53) — ( 54) 

N2 

8 

N2 

(54) 

X 

8 

ka 

(52) — (55 ) 

P 

8 

*o 

(52 ) — ( 54 ) 

YO 

12 

POm  and  QOm 

(7), (9), (12), (15), (14), (16) 

BJ (N+l ) 

21 

JN(ka) 

(53)— (55) 

BY  (N+l ) 

21 

YN(ka) 

(53 ) — ( 55 ) 

BJ1(N) 

26 

VNV 

(53 ) — ( 54 ) 

U1 

30 

l/(4H^2)(ka)) 

(53) 

U2 

31 

1/(8J  (ka)H(2) (ka)) 
o  o 

j,  m  ) 

(54) 

B 

39 

I  o 

Ni 

0 

(53)— (54) 

BJN 

41 

(ka) 

(53 ) — (55 ) 

U3 

42-45 

J,(N6  )  m 

<  J*  >'4  <ks) 

0 

(53)— (54) 

U 

46 

.N 

J 

(53) 

U1 

47 

U1 

(53) 

N2 

U2 

48 

U2  withoxit  j  J 

(54) 

N=N1 

BJY(N) 

51 

JN(ka)YN(ka) 

(54) 

B 

57 

J,  (N<h  )/N6 

(54) 

CO  1  c 
002C 
C03 
004 
C05 
G06 
C07 
008 
009  17 
010 
01127 
012 
013  11 
014 
015  12 
Cl  6 
017 
Cl  8 
019 
020 
021 
022 
02  J  18 
024 
02519 
026 
027 
C2820 
029 
03  Q 
031 
C32 
033  15 
034 
03  5  25 
03  6 
037 
C38 
039 
04  0 
041 
042 
04322 
C44 
04523 
C4624 
04  7 
C48 
049 
C5013 
051 
052 
05321 


LISTING  CF  THE  MAIN  PROGRAM  FOB  THE 
F003IEB  SEBIES  METHOD  QF  SCLOTICN 
COMPLEX  01, 02. U, 03, V 

COMMON  Y  0  (  3  3)  ,  PI2  ,  PI4  ,  PI7  ,  B  J  { 1  00 )  ,.EY{100) 
DIMENSION  BJ1 (10000) ,BJY (10000) 

OPEN  (0NIT  =  20,FILE=' MA0TZ3.DAT') 

CPEN  10NIT=21,FILE=' MA0TZ4.DAT') 

R  EAD  (20,  17)  N1,N2,X,P 
FOB  MAT  (I3,I5,2E14.7) 

WRITE  (21, 27)  N1,N2,X,P 

FOR  NAT ( *  N1  = ' ,13  , ' ,  N2=',I5,»,  X=*,B14.7,' 
READ  (20,11)  (YO  (I)  ,1=1,33) 

FORMAT  (5E14.7) 

WRITE  (21,12)  (YO  (I)  ,1=1,33) 

FORMAT  { *  YO  •/ ( IX,  5Z  1  4. 7)  ) 

?I=3.  141593 

PI2=2./PI 

EI4=FI/4. 

PI7=, 75*PI 
N 1  M=  N 1  - 1 
CALL  EES ( N1 M  ,X) 

WRITE  (21,18)  (BJ(I)  ,1=1, N1) 

FORM  AT  (  *  BJ  V  (1X,5E14.7)  ) 

WRITE  (21,19)  (BY  (I)  ,1=1,  N1) 

FORM  AT  (  *  BY  V(1I,5E14.7)  ) 

CALL  EESJ1  (N2,P, EJ  1) 

WRIT E  (2 1 , 20)  BJ1  (1)  ,BJ1  (N2) 

FORMAT  ('  BJ1  (1)  =  «,E14.7,»,  EJ  1  (N2)  ='  ,  E14." 
N=0 

01=. 25/C MPLI  (BJ  ( 1)  ,  -E  Y  ( 1)  ) 

U2  =  .  5/BJ  ( 1)  *01 
WRITE  (21,15) 

FORMAT  ( *  N’,  14X,'U1*,26X,,02*) 

WRITE  (21,25)  N,01,U2 
FORMAT (  1X,I3,4E14.7) 

0=  (1.  ,0.) 

CO  13  N=  1 , N  IN 
FN=N*P 
B=B J  1  (N) /FN 
NP=  N  +  1 
EJN=EJ (NP) 

IF  (IBS  (BJN)-I.  E-10)  22,22,  23 
03=B/B Y  (NP)  *(0.,1.) 

GO  TO  24 

0  3=B/CHPLX(EJN,-BY  (NP) ) 

0=0*  (0.  ,  1.) 

0 1=0*03  +  0  1 
02=B/EJN*03+02 
WRITE  (21, 25)  N,U1 ,02 
CCNTINOE 

CALL  BES J Y(N1,N2,X,EJY) 

WRITE  (21,21)  BJY (Nl) , BJY  (N2) 

FORM  AT  ( ’  BJY(N1)=»,E14.7,»,  BJY  (N2)  =  *  ,514, 


C54 

C55 

C56 

057 

osa 

059 
06C  14 
061 
062  16 
C63 
064 
065 
066 
067  26 
068 
C69 


£•=0. 

CO  14  N=N1,N2 
FN=N*P 
E=BJ  1  00  /FH 
S=S*E*3/BJI  (N) 

EJ1  (N)  =S 
CONTINUE 

WHITE  (21,  16)  (BJ1  (I)  ,I=N1,N2,50) 
FCHBAT  {  '  BJ1  V(1X,6E11.4)  ) 

02=S*  (0.  ,  1.  )  *02 
S=X*COS  (P) 

V=-2./P*01/02*CaPIX (COS (S) ,-SIN (S) ) 
SB  IT  E  (2  1 , 26)  V 
FCHHAT  (*  7=',2E14.7) 

STOP 

END 


C  LISTING  OP  THE  INPOT  DATA  FILE  HA0TZ3.DAT 
20 10 C 00  0-  1570796E  +  01  0. 8726646 E-0 1 

■0. 3072582E*04  0.736e758E+04-0.6065 1COE*03  0-  17 1 0234E+02-0. 
0.  1600171E-02-0.5961C89E-05  0- 9545773E-08  0. 41 631  50E  +  05  0. 
0.  1C0O0COE*O1-O.6O24727E*O4  0.  1 6  135  12E* 04- C. 75322 10E+02  0. 
-0.  12756C2E-01  G.5832787E-04-0. 1 1 C7698E-06  C . 30729 46E+05  0. 
0.  1000CC0E*01  0.95S9999E+00-0. 1C97659E-02  0. 246 1  455E-04  0. 
0.  18298 93 E— 02-0. 3191328E -04-0.  1562498E-01  0.  1427C 79E-0 3-0. 
0- 4667498 E-0 1-0- 1996720E-03  C.7217495E-C5 


2  27 10  0  1 E  +  00 
34202 1 1E+03 
1402590E+01 
2886431E+03 
1 00000  0  E  +  O  1 
5937434E-C5 
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C  LISTING  CF  THE  CCTPUT  DATA  FILE  KA0TZ4.DAT 

C 

N1  =  20.  N2= 1 0  COO,  X=  0.  1570796E+01 ,  P=  0.  87 26646E-01 
YO 

-0.  307 2582  E  +  O  4  0.  73 6875 8E *04-0. 60  85 1 00E+03  0. 171 0234E+02-0. 227100 1E+00 
0.  1600171  E-02-0.5961089E-05  0. 9545773E-08  0. 4 1 63 1 50E+05  0. 34202 1 1E+03 
0.  1 000 OCOE +  01-0.  60 2 4 72 7E *04  0.  16 135 1 2E+ 04- C.75322  1 0E  +  02  0.  1 402590E+01 
-0.  12756C2E-01  0. 5832787E-04-0. 11 C7698E-06  0. 3072946E*05  0. 288643  1E  +  03 
0.  1000000E+01  0.99S9999E  +  00-0. 1C97659E-02  0.2461 455E-04  0. 1000000E+01 
0.  18298S3E-02-0.3191328E-0  4-0.  1562498E-01  0.  1 427079E-03-0. 59374 34E-05 
0.  4687498E-01-0. 1998720E-03  0. 73 17495E-05 
BJ 

0. 4720014E+00  0.5668241E+00  0. 2497C16E+00  0.6903585E-01  0. 1399603E-01 
0. 2245355E-02  0.2983472E-03  0.3385059E-04  0. 3352192E-05  0. 2945642E-06 
0.  2326610E-07  0.  1669C26E-08  0.  1C96726E-09  C. 66484 97E-1  1  0. 37408 1 4 E- 12 
0.  1963744E-13  0.9661458E-1  5  0. 44726 1 5E-16  0.  1955078E-1 7  0.  8094814E- 19 
EY 

0. 4100  Q35E  +  0  0-0.3662806E 4-00-0. 8  763665E+0 0-0.  1865369E+0 1-0. 62488 19E+01 
-0. 2995961 E* 02-0. 18448CCE+0 3-0. 1 3793 64E+ 04- 0.12 109 35E+05-0. 1 2 196 55 E+ 06 
-0.  1385  5 13E+  07-0.  175 1893E+08-0. 2439783E+09-0. 37 10 1 S6E+ 1 0-0. 6 1 1 67 62E+ 11 
-0.  10866  25E*13-0.20  69184E*14-0-4204441E*15-C.9079853E+16-0. 2076745E+18 
B J  1  { 1 ) =  C. 4359171E-01,  BJ1 (N2)=-0. 2690578 E- Cl 
N  01  02 

0  0.3  C18775E+00  0. 2622256E+00  0.3  1 97845E+00  0. 2777805E+00 

1  Q.7036024E+0C  0.  8839C02E  +  QC  C. 8676468E* CO-O. 76 24704E-0 1 

2  0. 55381 85E+0C  0.  1409590E+0 1  0. 1 166432E+C 1-0. 1 1 2488 1 E*01 

3  0. 28842 84E+ GO  0. 1399768E+0 1  C. 1236961E+01-0.3030578E+01 

4  0.2686049E+00  0.  1320966E  +  0  1  0.  1 243 1  70E+01-0. 58C3C8 1 E+O  1 

5  0.3C49000E+00  0.  132C 967E»0 1  0.  1243436E+  01-0.9346033 E+0 1 

6  G.3C49000E+00  0. 1323586 E*01  0.  1 243443E*01-0.  1358569E+02 

7  0. 30 4554 1 E  +  00  0.  1323586E  +  01  0.  1243443E* C1-0.1845940E+02 

8  0.304554 1E+00  0.  1323547E  +  0 1  0.  1243443E+C1-0.2390473E  +  02 

9  0. 3045579E+0G  0.  1323547E  +  C1  0. 1 243443E+ C 1-0. 29 85677 E+02 
10  0.  3  G45  579E+  00  0.  1323547E  +  0  1  0.  1 243443E+0 1-0. 3624775E+02 
1  1  0.  3045579E+00  0.  1323547E+01  0.  1243443E+  C 1-0. 43 00740 E+02 

12  0.3045579E+00  0.  1323547E+01  0.  124 3443E+  01-0- 50 06366 E+02 

13  0.3045579E+00  0.  1323547E-+ 0 1  0. 1243443E+ 01-0.5734352E+02 

14  C. 3 C45  579E+  00  0.  1323547E  +  0  1  0.  124 3443E+ 0 1-0. 6477395 E+02 
1  5  0.3045579E  +  00  0. 1323547E+01  0. 124 3443E+ 01 -0.7228288 E+02 

16  C.3045579E+00  0.  1323547E+01  0.  1243443E+  C1-0.7980009E+02 

17  0. 3045579E+00  0.  1323547E+C  1  0.  1243443E+  C 1 -0. 87258 19 E+02 

18  0.30 4557SE+00  0. 13 235 47E+0 1  0. 1 243443E+C 1-0.9459342E+02 

19  0.3045579E+00  0.  1323547E  +  01  0 . 1 243443E+ C 1-0. 1 0 17465 E+03 
BJY ( Nl) =-0.  15964  93 E-01,  BJY  {N2)  =  -0.  3183099E-04 

EJ  1 

-0.6917E+01-0.8147E+02-0.8840E+02-0.9237E+02-0.9453E+02-0.9555E+02 
-0.966  1E+02-0.9716E+C2-0.9763E+02-0.9807E+02-C.9831E+02-0. 9861E+02 
-0. 9881E+02-0.9898E+C2-0.99 16 E+0 2-0 . 9927E+02-C . 994 0E+02-0. 99 5 1 E+ 02 
-0. 9959 E+02 -0.9 96 9E+ 02- C.  99 75E+ 02- 0. 9983E+ 02- C . 9990E+02-0. 9994E+ 0 2 
-0.  1COOE+O3-Q.10COE+03-Q.  10 01 E+0 3-0.  1001E+03-0.  1  00  2E+03-0.  1  0  02E  +  0  3 
-0.  1002E *03-0.  ICC 3E  +  C3-C. 10C3E  +  03- 0. 10  03E  +  03-0. 1 0048+03 -0. 1 004 E  + 03 
-0. 1004E+03-Q. 10C4E+C3-0. 1004E+03-0. 1005E+03-0. 1005E+03-0. 10C5E+03 
-0.  1005E+03-C.  10C5E  +  C3-C.  1006E+ 03- 0. 1 0C6E+ 03-0.  1 006E+03-0.  1006 E+0 3 


0.  1006E+03-0.10C6E+C3-C.  10 06 E+03-0.  1 0061+03-0.  1 0075  +  0 3-0.  1 C07E+0 3 
0.  1007E  +  0  3-Q. 1QC7F  +  C3-C. 10C7E+G3-0. 10G7E+Q3-0.  1007E+03-0.  1007E  +  03 
0.  1007E+03-0. 1 OC7E+C3-0.  1 008 E+0 3-0. 1 008 E+0 3-0.  10Q8E+03-C.  1 C08E+03 
0.  1008E+0  3-0. 1CC8E  +  C3-C. 1 0 08E  + 0 3- C. 1 00  8E  + 03-0. 1 0 08E+03 -0. 1 0 08 E  +  0 3 
0.  1008  E+0  3-0.  1  C08E+03-0.  1008E+03-0.  1008E+03-Q.  1008E+03-C.  1  0C8E+33 
0- 1009 E+0 3-0. 1CC9E+03- 0. 10 C9E+03- 0. 1 009E+ 03-0. 1 00 9E+03-0. 1009 E+0 3 
0. 1CC9E  +  03- 0.1 0C9E+ 03-0.  1009E+03-0. 1 009E  +  0  3-C.  100  9E  +  03-0.  1G09E  +  0  3 
0.  1009E+03-0.  1 0 C9E  +  C3-  0.  1009E +03- 0. 1 009E+03-C. 1 009E+03-0.  1009 E+0 3 
0. 1CC9E+03-0. 10C9E+03-0. 10 C 9E+ 03-0. 1 009E+ 03-0. 1 009E+03-0. 1 009E+03 
0. 1009 E+0 3-0. 10C9E+03-0. 100 9E +03-0. 1 009E+03-0. 1 0 1 0E+03- C. 1 0 10E+03 
0. 1010E+03-0. 1010E+C3-0. 10 1 0E+ 03- 0. 1 0 1 CE+03-0. 1 01 CE+03-0. 1 0 10 E+0 3 
■0.  10  10E  +  03-C.1010E  +  C3-0. 10 1 0E  +  03-0. 1 0 1 0E+03-0. 1 0 1 0E  +  03-C . 1 0 10E  +  03 
■0.  10  10E+03-C.  10  10E  +  C3-0.  10  10E  +  03-0.  1 0  1 0E  +  03-0.  1  0  1  CE  +  03-0.  1 0  10E+03 
0.  1010E  +  03-0. 10 10E  +  03-0.  10 10 E  +  03-0. 1 0 1 0 E  +  0 3-0. 1010E  +  03-0. 1010E  +  0  3 
0.  10  10E+0  3-0.  10  10E  +  C3-0.  10 10E  +  03- 0. 1 01 CE  + 03-0. 1 0 1 0E  +  03-0. 1 01 0E+03 
0-10 10E+03-0.1010E+03-0. 10 10E+03-0. 1010E+03-C. 1 0 1 0E+0 3- C. 1 0 1 0E+0 3 
0.  10 10 E+0 3-0. 10  10 E  +  C 3-0. 10 10E  +  03-0. 1 01 OE+03-0. 1 01 0E  +  03-0. 1010E  +  03 
0.  1010E+0 3-0. 1010 E+C3-0. 10  10 E+0 3-0. 10 1 0E  +  0  3-0.  10 1 0E  +  03-C. 1 0 1 0E+ 03 
0. 1010E+03-0. 1010E+C3-0- 10 1 0E+03-0. 1 0 1 0E+03-0. 1 01 0E+03 -0. 1010E+03 
0.  1010E+03-0. 10 11E+ 03-0.  10  11  E+0  3-0.  101 1  E+03-0.  101  1E  +  03-0.  10  1 1E+0  3 
0. 10 11 E+03-0- 10 11E+C3-0. 1011E+03-0. 101 1E+03-C. 1 011E+03-0. 1011E+03 
0. 101 1 E+03-0. 101 1E+ 03-0. 10 11 E+03-0. 1011E  +  Q3-0. 1 0 1  1E  +  0 3-0. 1 0 1 1E+0  3 
0. 1011E+03-0.1011E+03-0. 10 1 1 E+03- 0. 1 0 1 1 E+03-0 . 1 01 1 E+ 03-0. 1 0 1 1 E+0 3 
0.  1011 E+03-0. 1 01 1S  +  03-0. 10 1 1 E+0  3-0. 1 0 1 1 E+0 3-0. 10 11E  +  03-C. 101 1E  +  0  3 
0. 10 11 E+0 3-0. 10 11 E+03- C. 10 1 1E+03-0. 1 0 1 1E+03- 0. 1 01 1E+03-0. 1 0 1 1E+0 3 
0.  1011E  +  03-0. 1011E+03 
=-0. 34435  C  0E-01-0.  1495431E  +  00 
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